FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
hecmw_ordering_qmd.f90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
2 ! Copyright (c) 2020 FrontISTR Commons
3 ! This software is released under the MIT License, see LICENSE.txt
4 !-------------------------------------------------------------------------------
5 
6 !----------------------------------------------------------------------
8 ! ordering using quotient graphs
9 !----------------------------------------------------------------------
11  use hecmw_util
12  implicit none
13 
14  private
15  public :: hecmw_ordering_genqmd
16 
17 contains
18 
19  !======================================================================!
21  !======================================================================!
22  subroutine hecmw_ordering_genqmd(Neqns,Nttbr,Xadj,Adj0,Perm,Invp)
23  implicit none
24  !------
25  integer(kind=kint), intent(in):: neqns
26  integer(kind=kint), intent(in):: nttbr
27  integer(kind=kint), intent(in):: adj0(:)
28  integer(kind=kint), intent(in):: xadj(:)
29  integer(kind=kint), intent(out):: perm(:)
30  integer(kind=kint), intent(out):: invp(:)
31  !------
32  integer(kind=kint), allocatable:: deg(:)
33  integer(kind=kint), allocatable:: marker(:)
34  integer(kind=kint), allocatable:: rchset(:)
35  integer(kind=kint), allocatable:: nbrhd(:)
36  integer(kind=kint), allocatable:: qsize(:)
37  integer(kind=kint), allocatable:: qlink(:)
38  integer(kind=kint), allocatable:: adjncy(:)
39  integer(kind=kint):: inode
40  integer(kind=kint):: ip
41  integer(kind=kint):: irch
42  integer(kind=kint):: j
43  integer(kind=kint):: mindeg
44  integer(kind=kint):: ndeg
45  integer(kind=kint):: nhdsze
46  integer(kind=kint):: node
47  integer(kind=kint):: np
48  integer(kind=kint):: num
49  integer(kind=kint):: nump1
50  integer(kind=kint):: nxnode
51  integer(kind=kint):: rchsze
52  integer(kind=kint):: search
53  integer(kind=kint):: thresh
54  integer(kind=kint):: ierror
55  logical:: found
56 
57  allocate (deg(neqns+1),stat=ierror)
58  if ( ierror/=0 ) stop "ALLOCATION ERROR, deg: SUB. genqmd"
59  allocate (marker(neqns+1),stat=ierror)
60  if ( ierror/=0 ) stop "ALLOCATION ERROR, marker: SUB. genqmd"
61  allocate (rchset(neqns+2),stat=ierror)
62  if ( ierror/=0 ) stop "ALLOCATION ERROR, rchset: SUB. genqmd"
63  allocate (nbrhd(neqns+1),stat=ierror)
64  if ( ierror/=0 ) stop "ALLOCATION ERROR, nbrhd: SUB. genqmd"
65  allocate (qsize(neqns+1),stat=ierror)
66  if ( ierror/=0 ) stop "ALLOCATION ERROR, qsize: SUB. genqmd"
67  allocate (qlink(neqns+1),stat=ierror)
68  if ( ierror/=0 ) stop "ALLOCATION ERROR, qlink: SUB. genqmd"
69  allocate (adjncy(2*nttbr),stat=ierror)
70  if ( ierror/=0 ) stop "ALLOCATION ERROR, adjncy: SUB. genqmd"
71 
72  mindeg = neqns
73  adjncy(1:xadj(neqns+1) - 1) = adj0(1:xadj(neqns+1) - 1)
74  do node = 1, neqns
75  perm(node) = node
76  invp(node) = node
77  marker(node) = 0
78  qsize(node) = 1
79  qlink(node) = 0
80  ndeg = xadj(node+1) - xadj(node)
81  deg(node) = ndeg
82  if ( ndeg<mindeg ) mindeg = ndeg
83  enddo
84 
85  num = 0
86  loop1: do
87  search = 1
88  thresh = mindeg
89  mindeg = neqns
90  loop2: do
91  nump1 = num + 1
92  if ( nump1>search ) search = nump1
93  found = .false.
94  do j = search, neqns
95  node = perm(j)
96  if ( marker(node)>=0 ) then
97  ndeg = deg(node)
98  if ( ndeg<=thresh ) then
99  found = .true.
100  exit
101  endif
102  if ( ndeg<mindeg ) mindeg = ndeg
103  endif
104  enddo
105  if (.not. found) cycle loop1
106 
107  search = j
108  marker(node) = 1
109  call qmdrch(node,xadj,adjncy,deg,marker,rchsze,rchset,nhdsze,nbrhd)
110  nxnode = node
111  do
112  num = num + 1
113  np = invp(nxnode)
114  ip = perm(num)
115  perm(np) = ip
116  invp(ip) = np
117  perm(num) = nxnode
118  invp(nxnode) = num
119  deg(nxnode) = -1
120  nxnode = qlink(nxnode)
121  if ( nxnode<=0 ) then
122  if ( rchsze>0 ) then
123  call qmdupd(xadj,adjncy,rchsze,rchset,deg,qsize,qlink,marker,rchset(rchsze+1:),nbrhd(nhdsze+1:))
124  marker(node) = 0
125  do irch = 1, rchsze
126  inode = rchset(irch)
127  if ( marker(inode)>=0 ) then
128  marker(inode) = 0
129  ndeg = deg(inode)
130  if ( ndeg<mindeg ) mindeg = ndeg
131  if ( ndeg<=thresh ) then
132  mindeg = thresh
133  thresh = ndeg
134  search = invp(inode)
135  endif
136  endif
137  enddo
138  if ( nhdsze>0 ) call qmdot(node,xadj,adjncy,marker,rchsze,rchset,nbrhd)
139  endif
140  if ( num>=neqns ) exit
141  cycle loop2
142  endif
143  enddo
144  exit
145  enddo loop2
146  exit
147  enddo loop1
148 
149  deallocate (deg)
150  deallocate (marker)
151  deallocate (rchset)
152  deallocate (nbrhd)
153  deallocate (qsize)
154  deallocate (qlink)
155  deallocate (adjncy)
156  end subroutine hecmw_ordering_genqmd
157 
158  !======================================================================!
160  !======================================================================!
161  subroutine qmdrch(Root,Xadj,Adjncy,Deg,Marker,Rchsze,Rchset,Nhdsze,Nbrhd)
162  implicit none
163  !------
164  integer(kind=kint), intent(in):: root
165  integer(kind=kint), intent(in):: adjncy(:)
166  integer(kind=kint), intent(in):: deg(:)
167  integer(kind=kint), intent(in):: xadj(:)
168  integer(kind=kint), intent(out):: nhdsze
169  integer(kind=kint), intent(out):: rchsze
170  integer(kind=kint), intent(inout):: marker(:)
171  integer(kind=kint), intent(out):: rchset(:)
172  integer(kind=kint), intent(out):: nbrhd(:)
173  !------
174  integer(kind=kint):: i
175  integer(kind=kint):: istrt
176  integer(kind=kint):: istop
177  integer(kind=kint):: j
178  integer(kind=kint):: jstrt
179  integer(kind=kint):: jstop
180  integer(kind=kint):: nabor
181  integer(kind=kint):: node
182 
183  nhdsze = 0
184  rchsze = 0
185  istrt = xadj(root)
186  istop = xadj(root+1) - 1
187  if ( istop<istrt ) return
188  do i = istrt, istop
189  nabor = adjncy(i)
190  if ( nabor==0 ) return
191  if ( marker(nabor)==0 ) then
192  if ( deg(nabor)<0 ) then
193  marker(nabor) = -1
194  nhdsze = nhdsze + 1
195  nbrhd(nhdsze) = nabor
196  loop1: do
197  jstrt = xadj(nabor)
198  jstop = xadj(nabor+1) - 1
199  do j = jstrt, jstop
200  node = adjncy(j)
201  nabor = -node
202  if ( node<0 ) cycle loop1
203  if ( node==0 ) exit
204  if ( marker(node)==0 ) then
205  rchsze = rchsze + 1
206  rchset(rchsze) = node
207  marker(node) = 1
208  endif
209  enddo
210  exit
211  enddo loop1
212  else
213  rchsze = rchsze + 1
214  rchset(rchsze) = nabor
215  marker(nabor) = 1
216  endif
217  endif
218  enddo
219  end subroutine qmdrch
220 
221  !======================================================================!
223  !======================================================================!
224  subroutine qmdupd(Xadj,Adjncy,Nlist,List,Deg,Qsize,Qlink,Marker,Rchset,Nbrhd)
225  implicit none
226  !------
227  integer(kind=kint), intent(in):: nlist
228  integer(kind=kint), intent(in):: adjncy(:)
229  integer(kind=kint), intent(in):: list(:)
230  integer(kind=kint), intent(in):: xadj(:)
231  integer(kind=kint), intent(inout):: deg(:)
232  integer(kind=kint), intent(inout):: marker(:)
233  integer(kind=kint), intent(out):: rchset(:)
234  integer(kind=kint), intent(out):: nbrhd(:)
235  integer(kind=kint), intent(inout):: qsize(:)
236  integer(kind=kint), intent(inout):: qlink(:)
237  !------
238  integer(kind=kint):: deg0
239  integer(kind=kint):: deg1
240  integer(kind=kint):: il
241  integer(kind=kint):: inhd
242  integer(kind=kint):: inode
243  integer(kind=kint):: irch
244  integer(kind=kint):: j
245  integer(kind=kint):: jstrt
246  integer(kind=kint):: jstop
247  integer(kind=kint):: mark
248  integer(kind=kint):: nabor
249  integer(kind=kint):: nhdsze
250  integer(kind=kint):: node
251  integer(kind=kint):: rchsze
252 
253  if ( nlist<=0 ) return
254  deg0 = 0
255  nhdsze = 0
256  do il = 1, nlist
257  node = list(il)
258  deg0 = deg0 + qsize(node)
259  jstrt = xadj(node)
260  jstop = xadj(node+1) - 1
261  do j = jstrt, jstop
262  nabor = adjncy(j)
263  if ( marker(nabor)==0 .and. deg(nabor)<0 ) then
264  marker(nabor) = -1
265  nhdsze = nhdsze + 1
266  nbrhd(nhdsze) = nabor
267  endif
268  enddo
269  enddo
270 
271  if ( nhdsze>0 ) call qmdmrg(xadj,adjncy,deg,qsize,qlink,marker,deg0,nhdsze,nbrhd,rchset,nbrhd(nhdsze+1:))
272  do il = 1, nlist
273  node = list(il)
274  mark = marker(node)
275  if ( mark<=1 .and. mark>=0 ) then
276  call qmdrch(node,xadj,adjncy,deg,marker,rchsze,rchset,nhdsze,nbrhd)
277  deg1 = deg0
278  if ( rchsze>0 ) then
279  do irch = 1, rchsze
280  inode = rchset(irch)
281  deg1 = deg1 + qsize(inode)
282  marker(inode) = 0
283  enddo
284  endif
285  deg(node) = deg1 - 1
286  if ( nhdsze>0 ) then
287  do inhd = 1, nhdsze
288  inode = nbrhd(inhd)
289  marker(inode) = 0
290  enddo
291  endif
292  endif
293  enddo
294  end subroutine qmdupd
295 
296  !======================================================================!
298  !======================================================================!
299  subroutine qmdmrg(Xadj,Adjncy,Deg,Qsize,Qlink,Marker,Deg0,Nhdsze,Nbrhd,Rchset,Ovrlp)
300  implicit none
301  !------
302  integer(kind=kint), intent(in):: deg0
303  integer(kind=kint), intent(in):: nhdsze
304  integer(kind=kint), intent(in):: adjncy(:)
305  integer(kind=kint), intent(in):: nbrhd(:)
306  integer(kind=kint), intent(in):: xadj(:)
307  integer(kind=kint), intent(inout):: deg(:)
308  integer(kind=kint), intent(inout):: qsize(:)
309  integer(kind=kint), intent(inout):: qlink(:)
310  integer(kind=kint), intent(inout):: marker(:)
311  integer(kind=kint), intent(out):: rchset(:)
312  integer(kind=kint), intent(out):: ovrlp(:)
313  !------
314  integer(kind=kint):: deg1
315  integer(kind=kint):: head
316  integer(kind=kint):: inhd
317  integer(kind=kint):: iov
318  integer(kind=kint):: irch
319  integer(kind=kint):: j
320  integer(kind=kint):: jstrt
321  integer(kind=kint):: jstop
322  integer(kind=kint):: link
323  integer(kind=kint):: lnode
324  integer(kind=kint):: mark
325  integer(kind=kint):: mrgsze
326  integer(kind=kint):: nabor
327  integer(kind=kint):: node
328  integer(kind=kint):: novrlp
329  integer(kind=kint):: rchsze
330  integer(kind=kint):: root
331 
332  if ( nhdsze<=0 ) return
333  do inhd = 1, nhdsze
334  root = nbrhd(inhd)
335  marker(root) = 0
336  enddo
337  do inhd = 1, nhdsze
338  root = nbrhd(inhd)
339  marker(root) = -1
340  rchsze = 0
341  novrlp = 0
342  deg1 = 0
343  loop1: do
344  jstrt = xadj(root)
345  jstop = xadj(root+1) - 1
346  do j = jstrt, jstop
347  nabor = adjncy(j)
348  root = -nabor
349  if ( nabor<0 ) cycle loop1
350  if ( nabor==0 ) exit
351  mark = marker(nabor)
352 
353  if ( mark>=0 ) then
354  if ( mark<=0 ) then
355  rchsze = rchsze + 1
356  rchset(rchsze) = nabor
357  deg1 = deg1 + qsize(nabor)
358  marker(nabor) = 1
359  elseif ( mark<=1 ) then
360  novrlp = novrlp + 1
361  ovrlp(novrlp) = nabor
362  marker(nabor) = 2
363  endif
364  endif
365  enddo
366  exit
367  enddo loop1
368  head = 0
369  mrgsze = 0
370  loop2: do iov = 1, novrlp
371  node = ovrlp(iov)
372  jstrt = xadj(node)
373  jstop = xadj(node+1) - 1
374  do j = jstrt, jstop
375  nabor = adjncy(j)
376  if ( marker(nabor)==0 ) then
377  marker(node) = 1
378  cycle loop2
379  endif
380  enddo
381  mrgsze = mrgsze + qsize(node)
382  marker(node) = -1
383  lnode = node
384  do
385  link = qlink(lnode)
386  if ( link<=0 ) then
387  qlink(lnode) = head
388  head = node
389  exit
390  else
391  lnode = link
392  endif
393  enddo
394  enddo loop2
395  if ( head>0 ) then
396  qsize(head) = mrgsze
397  deg(head) = deg0 + deg1 - 1
398  marker(head) = 2
399  endif
400  root = nbrhd(inhd)
401  marker(root) = 0
402  if ( rchsze>0 ) then
403  do irch = 1, rchsze
404  node = rchset(irch)
405  marker(node) = 0
406  enddo
407  endif
408  enddo
409  end subroutine qmdmrg
410 
411  !======================================================================!
413  !======================================================================!
414  subroutine qmdot(Root,Xadj,Adjncy,Marker,Rchsze,Rchset,Nbrhd)
415  implicit none
416  !------
417  integer(kind=kint), intent(in):: rchsze
418  integer(kind=kint), intent(in):: root
419  integer(kind=kint), intent(in):: marker(:)
420  integer(kind=kint), intent(in):: rchset(:)
421  integer(kind=kint), intent(in):: nbrhd(:)
422  integer(kind=kint), intent(in):: xadj(:)
423  integer(kind=kint), intent(inout):: adjncy(:)
424  !------
425  integer(kind=kint):: inhd
426  integer(kind=kint):: irch
427  integer(kind=kint):: j
428  integer(kind=kint):: jstrt
429  integer(kind=kint):: jstop
430  integer(kind=kint):: link
431  integer(kind=kint):: nabor
432  integer(kind=kint):: node
433 
434  irch = 0
435  inhd = 0
436  node = root
437  loop1: do
438  jstrt = xadj(node)
439  jstop = xadj(node+1) - 2
440  if ( jstop>=jstrt ) then
441  do j = jstrt, jstop
442  irch = irch + 1
443  adjncy(j) = rchset(irch)
444  if ( irch>=rchsze ) exit loop1
445  enddo
446  endif
447  link = adjncy(jstop+1)
448  node = -link
449  if ( link>=0 ) then
450  inhd = inhd + 1
451  node = nbrhd(inhd)
452  adjncy(jstop+1) = -node
453  endif
454  enddo loop1
455 
456  adjncy(j+1) = 0
457  do irch = 1, rchsze
458  node = rchset(irch)
459  if ( marker(node)>=0 ) then
460  jstrt = xadj(node)
461  jstop = xadj(node+1) - 1
462  do j = jstrt, jstop
463  nabor = adjncy(j)
464  if ( marker(nabor)<0 ) then
465  adjncy(j) = root
466  exit
467  endif
468  enddo
469  endif
470  enddo
471  end subroutine qmdot
472 
473 end module hecmw_ordering_qmd
hecmw_ordering_qmd
HECMW_ORDERING_QMD is a program for the minimum degree.
Definition: hecmw_ordering_qmd.f90:10
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
hecmw_ordering_qmd::hecmw_ordering_genqmd
subroutine, public hecmw_ordering_genqmd(Neqns, Nttbr, Xadj, Adj0, Perm, Invp)
hecmw_ordering_GENQMD
Definition: hecmw_ordering_qmd.f90:23