QR_MUMPS
dqrm_do_subtree.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 
46 subroutine dqrm_do_subtree(qrm_mat, fnum, flops)
47 
48  use qrm_mem_mod
49  use dqrm_spmat_mod
50  use dqrm_fdata_mod
51  use qrm_adata_mod
52  use dqrm_utils_mod
53  use qrm_common_mod
54  use qrm_sort_mod
56  implicit none
57 
58  type(dqrm_spmat_type), target :: qrm_mat
59  integer :: fnum
60  real(kind(1.d0)) :: flops
61 
62  type(dqrm_front_type), pointer :: front, cfront
63  type(qrm_adata_type), pointer :: adata
64  type(dqrm_fdata_type), pointer :: fdata
65 
66  integer, allocatable :: iwork(:)
67  real(kind(1.d0)), allocatable :: rwork(:)
68  integer :: node, p, c, i, j, cnt
69 
70  ! error management
71  integer :: err_act
72  character(len=*), parameter :: name='qrm_do_subtree'
73 
74  call qrm_err_act_save(err_act)
75 
76  ! simplify
77  adata => qrm_mat%adata
78  fdata => qrm_mat%fdata
79 
80  node = fnum
81  front => fdata%front_list(node)
82  ! write(*,*)'--->',fnum,20*front%n*front%nb
83  call qrm_aalloc(rwork, adata%rc(fnum)*qrm_mat%icntl(qrm_nb_))
84 
85  cnt = 0
86  call qrm_aalloc(iwork, qrm_mat%n)
87  __qrm_check_ret(name,'qrm_aalloc',9999)
88 
89  subtree: do
90 
91  front => fdata%front_list(node)
92 
93  if (front%status .eq. qrm_ready_) then
94  ! the front is ready to be assembled and factorized
95 
96  ! assemble
97  call dqrm_init_front(qrm_mat, node, .false., work=iwork)
98  __qrm_check_ret(name,'qrm_init_front',9999)
99 
100  ! clean
101  do p = adata%childptr(node), adata%childptr(node+1)-1
102  c = adata%child(p)
103  cfront => fdata%front_list(c)
104  call dqrm_clean_front(qrm_mat, cfront%num)
105  __qrm_check_ret(name,'qrm_clean_front',9999)
106  cnt = cnt+1
107  cfront%status = qrm_done_
108  end do
109 
110  ! factorize
111  call sequential_fct(front, flops)
112  __qrm_check_ret(name,'sequential_fct',9999)
113 
114  front%status = qrm_factorized_
115 
116  if(node .eq. fnum) exit subtree ! we're done
117 
118  node = adata%parent(node)
119  if(node .eq. 0) exit subtree
120 
121  else
122  ! go down the subtree
123  node = adata%child(adata%childptr(node+1)+front%status)
124  cycle subtree
125  end if
126 
127 
128  end do subtree
129 
130  !$ call omp_set_lock(qrm_mat%fdata%lock)
131  qrm_mat%fdata%done = qrm_mat%fdata%done+cnt
132  !$ call omp_unset_lock(qrm_mat%fdata%lock)
133 
134 
135  call qrm_adealloc(rwork)
136  call qrm_adealloc(iwork)
137  __qrm_check_ret(name,'qrm_adealloc',9999)
138 
139  call qrm_err_act_restore(err_act)
140  return
141 
142 9999 continue ! error management
143  call qrm_err_act_restore(err_act)
144  if(err_act .eq. qrm_abort_) then
145  call qrm_err_check()
146  end if
147 
148  return
149 
150 contains
151 
152 
153 
154 
155 
156 
157 
158 
159 !==========================================================================================
160 !==========================================================================================
161  subroutine sequential_fct(front, flops)
162  ! this routine does the sequential factorization of a front
163 
164  implicit none
165 
166  type(dqrm_front_type) :: front
167  real(kind(1.d0)) :: flops
168 
169  integer :: ib, i, j, m, n, k, thn, info
170 
171  ! error management
172  integer :: err_act
173  character(len=*), parameter :: name='sequential_fct'
174 
175  call qrm_err_act_save(err_act)
176 
177 
178  call qrm_get(qrm_mat, 'qrm_ib', ib)
179 
180  call qrm_arealloc(rwork, front%n*front%nb)
181  __qrm_check_ret(name,'qrm_aalloc',9999)
182 
183  do i=1, front%ne, ib
184  n = min(ib, front%n-i+1)
185  m = max(front%stair(min(i+n-1,front%ne))- i+1,0) ! the height of the panel
186  k = min(m,n)
187  j = i+n
188 
189  call dgeqrf(m, n, front%front(i,i), front%m, front%tau(i), rwork(1), &
190  & front%nb**2, info)
191  call dlarft( 'f', 'c', m, k, front%front(i,i), front%m, front%tau(i), &
192  & front%t(1,i), front%ib )
193 
194  flops = flops+qrm_count_flops(m, n, n, 'panel')
195 
196  n = front%n-i-n+1
197  if(n.gt.0) call dlarfb('l', 't', 'f', 'c', m, n, k, front%front(i,i), &
198  & front%m, front%t(1,i), front%ib, &
199  & front%front(i,j), front%m, rwork(1), front%n)
200  flops = flops+qrm_count_flops(m, n, k, 'update')
201 
202  ! do j=i+n,front%n, front%nb
203  ! n = min(front%nb,front%n-j+1)
204  ! call dlarfb('l', 't', 'f', 'c', m, n, k, front%front(i,i), &
205  ! & front%m, front%t(1,1,(i+ib-1)/ib), front%nb, &
206  ! & front%front(i,j), front%m, rwork(1), front%nb)
207  ! flops = flops+qrm_count_flops(m, n, k, 'update')
208  ! end do
209 
210  end do
211 
212  __qrm_check_ret(name,'qrm_adealloc',9999)
213 
214  call qrm_err_act_restore(err_act)
215  return
216 
217 9999 continue ! error management
218  call qrm_err_act_restore(err_act)
219  if(err_act .eq. qrm_abort_) then
220  call qrm_err_check()
221  end if
222 
223  return
224 
225  end subroutine sequential_fct
226 
227 
228 end subroutine dqrm_do_subtree
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.
This module contains routines for sorting.
This module contains the interfaces of all non-typed routines.
subroutine dqrm_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...
Generif interface for the ::dqrm_pgeti, ::dqrm_pgetr and.
This module contains the definition of the analysis data type.
integer, parameter qrm_done_
The data structure meant to store all the results of the factorization phase.
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine dqrm_init_front(qrm_mat, fnum, par, work)
This routine initializes a front.
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
integer, parameter qrm_factorized_
subroutine dqrm_clean_front(qrm_mat, fnum)
This routine performs the cleaning of a front.
Generic interface for the ::qrm_count_realflops ::qrm_count_pureflops.
This type defines the data structure used to store a matrix.
This type defines a data structure containing all the data related to a front.
This module contains the definition of all the data related to the factorization phase.
Generic interface for the qrm_arealloc_i qrm_arealloc_s qrm_arealloc_d qrm_arealloc_c qrm_arealloc_z...
integer, parameter qrm_ready_
subroutine sequential_fct(front, flops)
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38