QR_MUMPS
qrm_activate_front.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 
38 
52 subroutine _qrm_activate_front(qrm_mat, fnum, flops)
53 
54  use qrm_mem_mod
55  use _qrm_spmat_mod
56  use _qrm_fdata_mod
57  use _qrm_utils_mod
58  use qrm_common_mod
59  use qrm_sort_mod
61  implicit none
62 
63  type(_qrm_spmat_type), target :: qrm_mat
64  integer :: fnum
65  real(kind(1.d0)) :: flops
66 
67  integer :: n_rows_orig, i, j, row, col, c, child, roff, crows, f
68  integer :: m, n, np, rsize, hsize, mr, mc, mn, k, p, ne
69  integer :: frow, fcol, nb, b, father, thn
70  type(_qrm_front_type), pointer :: front
71  type(qrm_adata_type), pointer :: adata
72 
73  ! error management
74  integer :: err_act
75  character(len=*), parameter :: name='qrm_activate_front'
76 
77  call qrm_err_act_save(err_act)
78 
79 
80  thn = 0
81  !$ thn=omp_get_thread_num()
82 
83  front => null()
84 
85  ! to make things easier
86  front => qrm_mat%fdata%front_list(fnum)
87  father = qrm_mat%adata%parent(fnum)
88  adata => qrm_mat%adata
89 
90  ! sweep over the children to check whether there are small subtrees
91  ! that have to be treated
92  do p = adata%childptr(fnum), adata%childptr(fnum+1)-1
93  c = adata%child(p)
94  if(adata%small(c) .eq. 1) call _qrm_do_subtree(qrm_mat, c, flops)
95  __qrm_check_ret(name,'qrm_do_subtree',9999)
96  end do
97 
98  call _qrm_init_front(qrm_mat, fnum, .true.)
99  __qrm_check_ret(name,'qrm_init_front',9999)
100 
101  front%owner = thn
102 
103  ! clean
104  do p = adata%childptr(fnum), adata%childptr(fnum+1)-1
105  c = adata%child(p)
106  if(adata%small(c) .eq. 1) then
107  call _qrm_clean_front(qrm_mat, c)
108  __qrm_check_ret(name,'qrm_clean_front',9999)
109 
110  !$ call omp_set_lock(qrm_mat%fdata%lock)
111  qrm_mat%fdata%done = qrm_mat%fdata%done+1
112  !$ call omp_unset_lock(qrm_mat%fdata%lock)
113 
114  end if
115 
116  end do
117 
118  call qrm_err_act_restore(err_act)
119  return
120 
121 9999 continue ! error management
122  call qrm_err_act_restore(err_act)
123  if(err_act .eq. qrm_abort_) then
124  call qrm_err_check()
125  end if
126 
127  return
128 
129 end subroutine _qrm_activate_front
130 
131 
132 
133 
134 
135 
137 
147 subroutine _qrm_clean_front(qrm_mat, fnum)
149  use qrm_mem_mod
150  use _qrm_spmat_mod
151  use _qrm_fdata_mod
153  implicit none
154 
155  type(_qrm_spmat_type), target :: qrm_mat
156  integer :: fnum
157 
158 
159  type(_qrm_front_type), pointer :: front
160  integer :: i, j
161 
162  ! error management
163  integer :: err_act
164  character(len=*), parameter :: name='qrm_clean_front'
165 
166  call qrm_err_act_save(err_act)
167 
168  front => qrm_mat%fdata%front_list(fnum)
169  ! if(front%num.eq.3) then
170  ! open(4,file="ffr.m",action="write")
171  ! ! open(5,file="ffi.m",action="write")
172  ! write(4,'("Ffr=[")',advance='no')
173  ! ! write(5,'("Ffi=[")',advance='no')
174  ! do i=1, front%m
175  ! do j=1, front%n
176  ! write(4,'(f20.15,x)',advance='no')real(front%front(i,j))
177  ! ! write(5,'(f20.15,x)',advance='no')aimag(front%front(i,j))
178  ! end do
179  ! write(4,'(" ")')
180  ! ! write(5,'(" ")')
181  ! end do
182  ! write(4,'("];")')
183  ! ! write(5,'("];")')
184  ! close(4)
185  ! ! close(5)
186  ! end if
187 
188  call qrm_adealloc(front%rowmap)
189  call qrm_adealloc(front%ptable)
190  call qrm_adealloc(front%colmap)
191  __qrm_check_ret(name,'qrm_adealloc',9999)
192 
193  if(qrm_mat%icntl(5) .gt. 0) then
194  ! h has to be stored
195  call _qrm_store_h(front)
196  __qrm_check_ret(name,'qrm_store_h',9999)
197  end if
198  call qrm_adealloc(front%t)
199  __qrm_check_ret(name,'qrm_adealloc',9999)
200 
201  call _qrm_store_r(front)
202  __qrm_check_ret(name,'qrm_store_r',9999)
203 
204  call qrm_adealloc(front%front)
205  __qrm_check_ret(name,'qrm_adealloc',9999)
206 
207  call qrm_err_act_restore(err_act)
208  return
209 
210 9999 continue ! error management
211  call qrm_err_act_restore(err_act)
212  if(err_act .eq. qrm_abort_) then
213  call qrm_err_check()
214  end if
215 
216  return
217 
218 end subroutine _qrm_clean_front
219 
220 
221 #if defined (RFPF)
222 
223 ! !==========================================================================================
224 ! !==========================================================================================
225 ! subroutine _qrm_store_h(front)
226 ! ! this subroutine is meant to store the Householder vectors
227 ! ! computed during the facto fron front%front to front%h
228 ! use qrm_mem_mod
229 ! use _qrm_fdata_mod
230 ! use _qrm_rfpf_mod
231 ! implicit none
232 
233 ! type(_qrm_front_type) :: front
234 
235 ! integer :: i, j, cnt, p, c, tsize, cs, hsize
236 
237 ! ! error management
238 ! integer :: err_act
239 ! character(len=*), parameter :: name='qrm_store_h'
240 
241 ! call qrm_err_act_save(err_act)
242 
243 ! hsize=0
244 ! j = 0
245 ! do
246 ! cs = min(front%nb, front%ne-j)
247 ! if(cs .le. 0) exit
248 ! j = j+cs
249 ! hsize = hsize+ cs*(cs+1)/2 + (max(front%stair(j),j) -j)*cs
250 ! end do
251 
252 ! call qrm_aalloc(front%h, hsize)
253 ! __QRM_CHECK_RET(name,'qrm_aalloc',9999)
254 ! cnt=1
255 ! do c=1, front%ne, front%nb
256 ! tsize = min(front%nb, front%ne-c+1)
257 ! ! store V1 (the triangular part)
258 ! call _qrm_to_rfpf('l', 'u', tsize, front%front(c,c), front%m, front%h(cnt))
259 ! ! store V2 (the rectangular part)
260 ! cnt = cnt+(tsize*(tsize+1))/2
261 ! do j = c, c+tsize-1
262 ! do i = c+tsize, front%stair(c+tsize-1)
263 ! front%h(cnt) = front%front(i,j)
264 ! cnt = cnt+1
265 ! end do
266 ! end do
267 ! end do
268 
269 ! call qrm_err_act_restore(err_act)
270 ! return
271 
272 ! 9999 continue ! error management
273 ! call qrm_err_act_restore(err_act)
274 ! if(err_act .eq. qrm_abort_) then
275 ! call qrm_err_check()
276 ! end if
277 
278 ! return
279 ! end subroutine _qrm_store_h
280 
281 
282 
283 ! !==========================================================================================
284 ! !==========================================================================================
285 ! subroutine _qrm_store_r(front)
286 ! ! this subroutine is meant to store the R factor
287 ! ! computed during the facto fron front%front to front%r
288 ! use qrm_mem_mod
289 ! use _qrm_fdata_mod
290 ! use _qrm_rfpf_mod
291 ! implicit none
292 
293 ! type(_qrm_front_type) :: front
294 
295 ! integer :: i, j, cnt, p, tsize, r, rsize
296 
297 ! ! error management
298 ! integer :: err_act
299 ! character(len=*), parameter :: name='qrm_store_r'
300 
301 ! call qrm_err_act_save(err_act)
302 
303 ! rsize = (front%npiv*(front%npiv+1)/2 + front%npiv*(front%n-front%npiv))
304 ! call qrm_aalloc(front%r, rsize)
305 ! __QRM_CHECK_RET(name,'qrm_aalloc',9999)
306 
307 ! cnt=1
308 
309 ! do r=1, front%npiv, front%nb
310 ! tsize = min(front%nb, front%npiv-r+1)
311 ! call _qrm_to_rfpf('u', 'n', tsize, front%front(r,r), front%m, front%r(cnt))
312 ! cnt = cnt+(tsize*(tsize+1))/2
313 ! do j=r+tsize, front%n
314 ! do i= r, r+tsize - 1
315 ! front%r(cnt) = front%front(i,j)
316 ! cnt = cnt+1
317 ! end do
318 ! end do
319 ! end do
320 
321 ! call qrm_err_act_restore(err_act)
322 ! return
323 
324 ! 9999 continue ! error management
325 ! call qrm_err_act_restore(err_act)
326 ! if(err_act .eq. qrm_abort_) then
327 ! call qrm_err_check()
328 ! end if
329 
330 ! return
331 ! end subroutine _qrm_store_r
332 
333 
334 #else
335 !==========================================================================================
336 !==========================================================================================
337 subroutine _qrm_store_h(front)
338  ! this subroutine is meant to store the Householder vectors
339  ! computed during the facto fron front%front to front%h
340  use qrm_mem_mod
341  use _qrm_fdata_mod
342  use _qrm_rfpf_mod
343  implicit none
344 
345  type(_qrm_front_type) :: front
346 
347  integer :: i, j, cnt, p, c, cs, hsize, m, pk, jp, k
348 
349  ! error management
350  integer :: err_act
351  character(len=*), parameter :: name='qrm_store_h'
352 
353  call qrm_err_act_save(err_act)
354 
355  hsize=0
356  front%hsize=0
357 
358  cnt=1
359  do jp = 1, front%ne, front%nb
360  pk = min(front%nb, front%ne-jp+1)
361  if(pk .le. 0) exit
362 
363  do j = jp, jp+pk-1, front%ib
364  k = min(front%ib, jp+pk - j)
365  if(k .le. 0) exit
366  m = max(front%stair(j+k-1),j+k-1) - j+1
367  hsize = hsize+k*m
368  end do
369  end do
370 
371  call qrm_aalloc(front%h, hsize)
372  __qrm_check_ret(name,'qrm_aalloc',9999)
373 
374  cnt=1
375  j = 0
376  p = 1
377 
378  outer: do jp = 1, front%ne, front%nb
379  pk = min(front%nb, front%ne-jp+1)
380  if(pk .le. 0) exit outer
381 
382  inner: do j = jp, jp+pk-1, front%ib
383  k = min(front%ib, jp+pk - j)
384  if(k .le. 0) exit inner
385  m = max(front%stair(j+k-1),j+k-1) - j+1
386  front%hsize = front%hsize+k*(k+1)/2 + k*(m-k)
387  do c=1, k
388  front%h(cnt:cnt+c-1) = front%t(1:c, j+c-1)
389  front%h(cnt+c:cnt+m-1) = front%front(j+c:j+m-1,j+c-1)
390  cnt = cnt+m
391  end do
392  end do inner
393  p = p+1
394  end do outer
395 
396  call qrm_err_act_restore(err_act)
397  return
398 
399 9999 continue ! error management
400  call qrm_err_act_restore(err_act)
401  if(err_act .eq. qrm_abort_) then
402  call qrm_err_check()
403  end if
404 
405  return
406 end subroutine _qrm_store_h
407 
408 
409 
410 !==========================================================================================
411 !==========================================================================================
412 subroutine _qrm_store_r(front)
413  ! this subroutine is meant to store the R factor
414  ! computed during the facto fron front%front to front%r
415  use qrm_mem_mod
416  use _qrm_fdata_mod
417  use _qrm_rfpf_mod
418  implicit none
419 
420  type(_qrm_front_type) :: front
421 
422  integer :: i, j, cnt, rsize, n, c, jp, k, pk
423 
424  ! error management
425  integer :: err_act
426  character(len=*), parameter :: name='qrm_store_r'
427 
428  call qrm_err_act_save(err_act)
429 
430  rsize=0
431  front%rsize=0
432 
433  cnt=1
434  do jp = 1, front%npiv, front%nb
435  pk = min(front%nb, front%npiv-jp+1)
436  if(pk .le. 0) exit
437 
438  do j = jp, jp+pk-1, front%ib
439  k = min(front%ib, jp+pk - j)
440  if(k .le. 0) exit
441  n = front%n-j+1
442  rsize = rsize+k*n
443  end do
444  end do
445 
446  call qrm_aalloc(front%r, rsize)
447  __qrm_check_ret(name,'qrm_aalloc',9999)
448 
449  front%rsize = front%npiv*(front%npiv+1)/2 + &
450  & front%npiv*(front%n-front%npiv)
451 
452  cnt=1
453  outer: do jp = 1, front%npiv, front%nb
454  pk = min(front%nb, front%npiv-jp+1)
455  if(pk .le. 0) exit outer
456 
457  inner: do j = jp, jp+pk-1, front%ib
458  k = min(front%ib, jp+pk - j)
459  if(k .le. 0) exit inner
460  n = front%n-j+1
461 
462  do c=1, n
463  front%r(cnt:cnt+k-1) = front%front(j:j+k-1,j+c-1)
464  cnt = cnt+k
465  end do
466  end do inner
467  end do outer
468 
469  call qrm_err_act_restore(err_act)
470  return
471 
472 9999 continue ! error management
473  call qrm_err_act_restore(err_act)
474  if(err_act .eq. qrm_abort_) then
475  call qrm_err_check()
476  end if
477 
478  return
479 end subroutine _qrm_store_r
480 #endif
481 
482 
This module contains generic interfaces for a number of auxiliary tools.
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_init_front(qrm_mat, fnum, par, work)
This routine initializes a front.
This module contains routines for sorting.
This module contains the interfaces of all non-typed routines.
subroutine _qrm_do_subtree(qrm_mat, fnum, flops)
This subroutine does the sequential factorization of an entire subtree.
This module contains all the generic interfaces for the typed routines in the factorization phase...
subroutine _qrm_store_h(front)
subroutine _qrm_clean_front(qrm_mat, fnum)
This routine performs the cleaning of a front.
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...
This module contains the definition of all the data related to the factorization phase.
subroutine _qrm_store_r(front)
subroutine _qrm_activate_front(qrm_mat, fnum, flops)
This routine activates a front.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
This module contains an implementation of some operations on triangular/trapezoidal matrices stored i...
This type defines a data structure containing all the data related to a front.