35 #include "qrm_common.h" 60 real(kind(1.d0)) :: flops
66 integer,
allocatable :: iwork(:)
67 _qrm_data,
allocatable :: rwork(:)
68 integer :: node, p, c, i, j, cnt
72 character(len=*),
parameter :: name=
'qrm_do_subtree' 74 call qrm_err_act_save(err_act)
77 adata => qrm_mat%adata
78 fdata => qrm_mat%fdata
81 front => fdata%front_list(node)
87 __qrm_check_ret(name,
'qrm_aalloc',9999)
91 front => fdata%front_list(node)
98 __qrm_check_ret(name,
'qrm_init_front',9999)
101 do p = adata%childptr(node), adata%childptr(node+1)-1
103 cfront => fdata%front_list(c)
105 __qrm_check_ret(name,
'qrm_clean_front',9999)
112 __qrm_check_ret(name,
'sequential_fct',9999)
116 if(node .eq. fnum)
exit subtree
118 node = adata%parent(node)
119 if(node .eq. 0)
exit subtree
123 node = adata%child(adata%childptr(node+1)+front%status)
131 qrm_mat%fdata%done = qrm_mat%fdata%done+cnt
137 __qrm_check_ret(name,
'qrm_adealloc',9999)
139 call qrm_err_act_restore(err_act)
143 call qrm_err_act_restore(err_act)
144 if(err_act .eq. qrm_abort_)
then 167 real(kind(1.d0)) :: flops
169 integer :: ib, i, j, m, n, k, thn, info
173 character(len=*),
parameter :: name=
'sequential_fct' 175 call qrm_err_act_save(err_act)
178 call qrm_get(qrm_mat,
'qrm_ib', ib)
181 __qrm_check_ret(name,
'qrm_aalloc',9999)
184 n = min(ib, front%n-i+1)
185 m = max(front%stair(min(i+n-1,front%ne))- i+1,0)
189 call _xgeqrf(m, n, front%front(i,i), front%m, front%tau(i), rwork(1), &
191 call _xlarft(
'f',
'c', m, k, front%front(i,i), front%m, front%tau(i), &
192 & front%t(1,i), front%ib )
197 if(n.gt.0)
call _xlarfb(
'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)
212 __qrm_check_ret(name,
'qrm_adealloc',9999)
214 call qrm_err_act_restore(err_act)
218 call qrm_err_act_restore(err_act)
219 if(err_act .eq. qrm_abort_)
then 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.
integer, parameter qrm_done_
integer, parameter qrm_ready_
This module contains all the generic interfaces for the typed routines in the factorization phase...
This module contains the definition of the analysis data type.
integer, parameter qrm_factorized_
subroutine _qrm_clean_front(qrm_mat, fnum)
This routine performs the cleaning of 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.
This type defines the data structure used to store a matrix.
Generic interface for the ::qrm_count_realflops ::qrm_count_pureflops.
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...
Generic interface for the qrm_arealloc_i qrm_arealloc_s qrm_arealloc_d qrm_arealloc_c qrm_arealloc_z...
This module contains the definition of all the data related to the factorization phase.
Generif interface for the ::_qrm_pgeti, ::_qrm_pgetr and.
subroutine sequential_fct(front, flops)
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
This type defines a data structure containing all the data related to a front.