FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
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
17contains
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
473end module hecmw_ordering_qmd
HECMW_ORDERING_QMD is a program for the minimum degree.
subroutine, public hecmw_ordering_genqmd(neqns, nttbr, xadj, adj0, perm, invp)
hecmw_ordering_GENQMD
I/O and Utility.
Definition: hecmw_util_f.F90:7