FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_matrix_ordering_MC.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
5
7 use hecmw_util
8 implicit none
9
10 private
12
13contains
14
15 subroutine hecmw_matrix_ordering_mc(N, indexL, itemL, indexU, itemU, &
16 perm_cur, ncolor_in, ncolor_out, COLORindex, perm, iperm)
17 implicit none
18 integer(kind=kint), intent(in) :: n
19 integer(kind=kint), intent(in) :: indexl(0:), indexu(0:)
20 integer(kind=kint), intent(in) :: iteml(:), itemu(:)
21 integer(kind=kint), intent(in) :: perm_cur(:)
22 integer(kind=kint), intent(in) :: ncolor_in
23 integer(kind=kint), intent(out) :: ncolor_out
24 integer(kind=kint), intent(out) :: colorindex(0:)
25 integer(kind=kint), intent(out) :: perm(:), iperm(:)
26 integer(kind=kint), allocatable :: iwk(:)
27 integer(kind=kint) :: nn_color, cntall, cnt, color
28 integer(kind=kint) :: i, inode, j, jnode
29 allocate(iwk(n))
30 iwk = 0
31 nn_color = n / ncolor_in
32 cntall = 0
33 colorindex(0) = 0
34 do color=1,n
35 cnt = 0
36 do i=1,n
37 inode = perm_cur(i)
38 if (iwk(inode) > 0 .or. iwk(inode) == -1) cycle
39 ! if (iwk(inode) == 0)
40 iwk(inode) = color
41 cntall = cntall + 1
42 perm(cntall) = inode
43 cnt = cnt + 1
44 if (cnt == nn_color) exit
45 if (cntall == n) exit
46 ! mark all connected and uncolored nodes
47 do j = indexl(inode-1)+1, indexl(inode)
48 jnode = iteml(j)
49 if (iwk(jnode) == 0) iwk(jnode) = -1
50 end do
51 do j = indexu(inode-1)+1, indexu(inode)
52 jnode = itemu(j)
53 if (jnode > n) cycle
54 if (iwk(jnode) == 0) iwk(jnode) = -1
55 end do
56 end do
57 colorindex(color) = cntall
58 if (cntall == n) then
59 ncolor_out = color
60 exit
61 end if
62 ! unmark all marked nodes
63 do i=1,n
64 if (iwk(i) == -1) iwk(i) = 0
65 end do
66 end do
67 deallocate(iwk)
68 ! make iperm
69 do i=1,n
70 iperm(perm(i)) = i
71 end do
72 end subroutine hecmw_matrix_ordering_mc
73
I/O and Utility.
Definition: hecmw_util_f.F90:7
subroutine, public hecmw_matrix_ordering_mc(n, indexl, iteml, indexu, itemu, perm_cur, ncolor_in, ncolor_out, colorindex, perm, iperm)