36 #include "qrm_common.h"
75 subroutine qrm_amalg_tree(n, parent, rowcount, porder, nvar, min_var, fill_thresh)
83 integer :: parent(:), rowcount(:), porder(:), nvar(:)
84 real(kind(1.d0)) :: fill_thresh
86 integer,
allocatable :: snodes_ptr(:)
87 integer :: nsnodes,
i, j, k, f
88 integer,
allocatable :: snodes_map(:), snodes_parent(:), kids(:), &
89 & snodes_merged(:), snodes_nvar(:), snodes_rc(:)
90 real(kind(1.d0)) :: ratio
91 integer(kind=8),
allocatable :: snodes_fill(:)
92 integer(kind=8) :: totfill, msize, fill
93 integer :: pv, fpv, rc, frc, col, fcol, tcol,
s, ns
97 character(len=*),
parameter :: name=
'qrm_amalg_tree'
103 __qrm_check_ret(name,
'qrm_postorder',9999)
113 if(f .ne. 0) kids(f) = kids(f) + 1
122 if ( (parent(j) .ne.
i) .or. &
123 & (rowcount(j) .ne. rowcount(
i)+1) .or. &
124 & (kids(
i) .gt. 1))
then
126 snodes_ptr(nsnodes) = k
129 snodes_ptr(nsnodes+1)=n+1
134 call move_alloc(kids, snodes_map)
142 pv = porder(snodes_ptr(
i))
143 snodes_nvar(
i) = snodes_ptr(
i+1)-snodes_ptr(
i)
144 snodes_rc(
i) = rowcount(pv)
145 do j=snodes_ptr(
i), snodes_ptr(
i+1)-1
147 snodes_map(porder(j)) =
i
149 f = parent(porder(j-1))
153 snodes_parent(
i) = snodes_map(f)
158 allocate(snodes_fill(nsnodes))
160 call move_alloc(snodes_map, snodes_merged)
164 snodes_loop:
do s = nsnodes-1, 1, -1
170 if(snodes_merged(f) .eq. 0)
exit
194 fcol = snodes_nvar(f)
196 fill = col*(frc+col-rc)
197 totfill = snodes_fill(f) + fill
199 if(tcol .lt. min_var)
then
201 else if (fill .eq. 0)
then
202 merge = tcol.lt.min_var*50
204 msize = (tcol)*(tcol+1)/2 + tcol*(frc-fcol)
206 ratio =
real(totfill, kind(1.d0))/
real(msize,kind(1.d0))
208 merge = (tcol .lt. min_var*4 .and. ratio .lt. fill_thresh*16) .or.&
209 & (tcol .lt. min_var*12 .and. ratio .lt. fill_thresh*2) .or. &
210 & (ratio .lt. fill_thresh)
216 snodes_fill(
s) = totfill
217 snodes_nvar(
s) = tcol
220 snodes_rc(
s) = col+frc
232 pv = porder(snodes_ptr(
i))
233 if(snodes_merged(
i) .eq. 0)
then
235 rowcount(pv) = snodes_rc(
i)
236 nvar(pv) = snodes_nvar(
i)
237 snodes_ptr(ns) = snodes_ptr(ns-1)+snodes_nvar(
i)
238 if(snodes_ptr(ns).gt.n)
exit
243 call move_alloc(snodes_merged, snodes_map )
245 pv = porder(snodes_ptr(
i))
247 f = parent(porder(snodes_ptr(
i+1)-1))
248 do j=snodes_ptr(
i)+1, snodes_ptr(
i+1)-1
250 snodes_map(porder(j)) = pv
251 rowcount(porder(j)) = 0
252 parent(porder(j)) = -pv
257 parent(pv) = snodes_map(f)
272 if(err_act .eq. qrm_abort_)
then
This module contains all the error management routines and 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.
subroutine qrm_amalg_tree(n, parent, rowcount, porder, nvar, min_var, fill_thresh)
This subroutine performs amalgamation on the assembly tree.
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.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
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.
Generic interface for the ::qrm_amalg_tree routine.
Generic interface for the ::qrm_postorder routine.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.