稀疏矩阵算法最小度数算法八之最小度数算法子程序QMDMRG

本章节介绍 the minimum degree algorithm ,我个人将其翻译为最小度数算法。这里我们看最小度数算法的子程序QMDMRG。

QMDMRG

这个子程序根据之前提到的原理确定不可区分的节点。设 C1,C2,R1,R2,Y 如同引理 5.4.6。这个子程序假定 C1 和 R1 事先已经给定 (C1 就是当前新找到的需要消元的节点形成的 supernode , R1 是其相邻的未消元的节点)。注意 R1 中的节点的 MARKER 已经设为 1

输入变量 DEG0 是 R1 中的节点数目。

这里存在很多个 C2,存放在 (NHDSZE, NBRHD) ,其中每一个 NBRHD(i) 均为一个已经消元的 supernode 。

循环 DO 1400 对于每一个 C2 使用条件,它首先确定集合 R2 - R1 ,存放在 (RCHSZE, RCHSET) 。循环 DO 600 循环确定交集 ,存放在 (NOVRLP, OVRLP) 。对于交集中的每一个节点,通过 DO 1100 循环确定能否满足条件( ),如果满足条件,那么这个节点会包括在一个 supernode 中,并放置在 QLINK 向量中,同时计算这个新的 supernode 的大小,并更新其 DEG 值。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
C****************************************************************          1.
C**************************************************************** 2.
C********** QMDMRG ..... QUOT MIN DEG MERGE *********** 3.
C**************************************************************** 4.
C**************************************************************** 5.
C 6.
C PURPOSE - THIS ROUTINE MERGES INDISTINGUISHABLE NODES IN 7.
C THE MINIMUM DEGREE ORDERING ALGORITHM. 8.
C IT ALSO COMPUTES THE NEW DEGREES OF THESE 9.
C NEW SUPERNODES. 10.
C 11.
C INPUT PARAMETERS - 12.
C (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. 13.
C DEG0 - THE NUMBER OF NODES IN THE GIVEN SET. 14.
C (NHDSZE, NBRHD) - THE SET OF ELIMINATED SUPERNODES 15.
C ADJACENT TO SOME NODES IN THE SET. 16.
C 17.
C UPDATED PARAMETERS - 18.
C DEG - THE DEGREE VECTOR. 19.
C QSIZE - SIZE OF INDISTINGUISHABLE NODES. 20.
C QLINK - LINKED LIST FOR INDISTINGUISHABLE NODES. 21.
C MARKER - THE GIVEN SET IS GIVEN BY THOSE NODES WITH 22.
C MARKER VALUE SET TO 1. THOSE NODES WITH DEGREE 23.
C UPDATED WILL HAVE MARKER VALUE SET TO 2. 24.
C 25.
C WORKING PARAMETERS - 26.
C RCHSET - THE REACHABLE SET. 27.
C OVRLP - TEMP VECTOR TO STORE THE INTERSECTION OF TWO 28.
C REACHABLE SETS. 29.
C 30.
C**************************************************************** 31.
C 32.
SUBROUTINE QMDMRG ( XADJ, ADJNCY, DEG, QSIZE, QLINK, 33.
1 MARKER, DEG0, NHDSZE, NBRHD, RCHSET, 34.
1 OVRLP ) 35.
C 36.
C**************************************************************** 37.
C 38.
INTEGER ADJNCY(1), DEG(1), QSIZE(1), QLINK(1), 39.
1 MARKER(1), RCHSET(1), NBRHD(1), OVRLP(1) 40.
INTEGER XADJ(1), DEG0, DEG1, HEAD, INHD, IOV, IRCH, 41.
1 J, JSTRT, JSTOP, LINK, LNODE, MARK, MRGSZE, 42.
1 NABOR, NHDSZE, NODE, NOVRLP, RCHSZE, ROOT 43.
C 44.
C**************************************************************** 45.
C 46.
C ------------------ 47.
C INITIALIZATION ... 48.
C ------------------ 49.
IF ( NHDSZE .LE. 0 ) RETURN 50.
DO 100 INHD = 1, NHDSZE 51.
ROOT = NBRHD(INHD) 52.
MARKER(ROOT) = 0 53.
100 CONTINUE 54.
C ------------------------------------------------- 55.
C LOOP THROUGH EACH ELIMINATED SUPERNODE IN THE SET 56.
C (NHDSZE, NBRHD). 57.
C ------------------------------------------------- 58.
DO 1400 INHD = 1, NHDSZE 59.
ROOT = NBRHD(INHD) 60.
MARKER(ROOT) = - 1 61.
RCHSZE = 0 62.
NOVRLP = 0 63.
DEG1 = 0 64.
200 JSTRT = XADJ(ROOT) 65.
JSTOP = XADJ(ROOT+1) - 1 66.
C ---------------------------------------------- 67.
C DETERMINE THE REACHABLE SET AND ITS INTERSECT- 68.
C ION WITH THE INPUT REACHABLE SET. 69.
C ---------------------------------------------- 70.
DO 600 J = JSTRT, JSTOP 71.
NABOR = ADJNCY(J) 72.
ROOT = - NABOR 73.
IF (NABOR) 200, 700, 300 74.
C 75.
300 MARK = MARKER(NABOR) 76.
IF ( MARK ) 600, 400, 500 77.
400 RCHSZE = RCHSZE + 1 78.
RCHSET(RCHSZE) = NABOR 79.
DEG1 = DEG1 + QSIZE(NABOR) 80.
MARKER(NABOR) = 1 81.
GOTO 600 82.
500 IF ( MARK .GT. 1 ) GOTO 600 83.
NOVRLP = NOVRLP + 1 84.
OVRLP(NOVRLP) = NABOR 85.
MARKER(NABOR) = 2 86.
600 CONTINUE 87.
C -------------------------------------------- 88.
C FROM THE OVERLAPPED SET, DETERMINE THE NODES 89.
C THAT CAN BE MERGED TOGETHER. 90.
C -------------------------------------------- 91.
700 HEAD = 0 92.
MRGSZE = 0 93.
DO 1100 IOV = 1, NOVRLP 94.
NODE = OVRLP(IOV) 95.
JSTRT = XADJ(NODE) 96.
JSTOP = XADJ(NODE+1) - 1 97.
DO 800 J = JSTRT, JSTOP 98.
NABOR = ADJNCY(J) 99.
IF ( MARKER(NABOR) .NE. 0 ) GOTO 800 100.
MARKER(NODE) = 1 101.
GOTO 1100 102.
800 CONTINUE 103.
C ----------------------------------------- 104.
C NODE BELONGS TO THE NEW MERGED SUPERNODE. 105.
C UPDATE THE VECTORS QLINK AND QSIZE. 106.
C ----------------------------------------- 107.
MRGSZE = MRGSZE + QSIZE(NODE) 108.
MARKER(NODE) = - 1 109.
LNODE = NODE 110.
900 LINK = QLINK(LNODE) 111.
IF ( LINK .LE. 0 ) GOTO 1000 112.
LNODE = LINK 113.
GOTO 900 114.
1000 QLINK(LNODE) = HEAD 115.
HEAD = NODE 116.
1100 CONTINUE 117.
IF ( HEAD .LE. 0 ) GOTO 1200 118.
QSIZE(HEAD) = MRGSZE 119.
DEG(HEAD) = DEG0 + DEG1 - 1 120.
MARKER(HEAD) = 2 121.
C -------------------- 122.
C RESET MARKER VALUES. 123.
C -------------------- 124.
1200 ROOT = NBRHD(INHD) 125.
MARKER(ROOT) = 0 126.
IF ( RCHSZE .LE. 0 ) GOTO 1400 127.
DO 1300 IRCH = 1, RCHSZE 128.
NODE = RCHSET(IRCH) 129.
MARKER(NODE) = 0 130.
1300 CONTINUE 131.
1400 CONTINUE 132.
RETURN 133.
END 134.

首先进行初始化,将相邻的已经消元的 supernode 的 MARKER 值设置为 0(QMDUPD 子程序里一开始将这些节点的 MARKER 值设为 - 1 了,本来可以直接在 QMDUPD 子程序里将 MARKER 值设为0的(应该也不用设,本来MARKER值就是 0) )。

1
2
3
4
5
6
7
8
C        ------------------                                               47.
C INITIALIZATION ... 48.
C ------------------ 49.
IF ( NHDSZE .LE. 0 ) RETURN 50.
DO 100 INHD = 1, NHDSZE 51.
ROOT = NBRHD(INHD) 52.
MARKER(ROOT) = 0 53.
100 CONTINUE 54.

遍历 (NHDSZE, NBRHD) 中每一个已经消元的 supernode,设置为 ROOT ,将其 MARKER 值设为 -1 。JSTRT 和 JSTOP 是 ROOT 的相邻节点。

1
2
3
4
5
6
7
8
9
10
11
12
C        -------------------------------------------------                55.
C LOOP THROUGH EACH ELIMINATED SUPERNODE IN THE SET 56.
C (NHDSZE, NBRHD). 57.
C ------------------------------------------------- 58.
DO 1400 INHD = 1, NHDSZE 59.
ROOT = NBRHD(INHD) 60.
MARKER(ROOT) = - 1 61.
RCHSZE = 0 62.
NOVRLP = 0 63.
DEG1 = 0 64.
200 JSTRT = XADJ(ROOT) 65.
JSTOP = XADJ(ROOT+1) - 1 66.

DO 600 遍历 ROOT 的相邻节点,设为 NABOR 。

ROOT 设为 - NABOR 。如果 NABOR 小于0,说明此时NABOR 为链接,则 GO TO 200 ,继续遍历;如果 NABOR 等于0,则说明此时所有相邻节点,遍历结束,GO TO 700 ;如果 NABOR 大于0,继续往下运行,GO TO 300

将 NABOR 的 MARKER 值设为 MARK ,如果 MARK 小于0 (因为 NABOR 是没有消元的节点,此时只有一种情况,就是 NABOR 合并到与其不可区分的节点中 ),则GO TO 600,跳过这个节点;如果 MARK 等于0,说明这个节点不在 R1 中,即其在 R2 - R1 中,则继续往下运行,则GO TO 400;如果 MARK 大于0,说明其同样在 R1 中,即其在 中, 则GO TO 500

我们先看 400 的语句,如果 MARK 等于0,确定集合 R2 - R1 (在R2 但不在 R1 中,因为 R1 的 MARKER 值均为 1) ,存放在 (RCHSZE, RCHSET),将其 MARKER 值设为 1 , DEG1 添加这些节点的数目。

我们再看 500 的语句,如果 MARK 等于2(说明这个节点已经是别的 C2 找到的不可区分的 supernode,这种情况感觉不可能发生,因为其相邻节点满足 ,因此不可能在别的 中 ),跳过这个节点。如果 MARK 等于1,确定为交集 ,存放在 (NOVRLP, OVRLP) ,将其 MARKER 值设为 2。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
C           ----------------------------------------------                67.
C DETERMINE THE REACHABLE SET AND ITS INTERSECT- 68.
C ION WITH THE INPUT REACHABLE SET. 69.
C ---------------------------------------------- 70.
DO 600 J = JSTRT, JSTOP 71.
NABOR = ADJNCY(J) 72.
ROOT = - NABOR 73.
IF (NABOR) 200, 700, 300 74.
C 75.
300 MARK = MARKER(NABOR) 76.
IF ( MARK ) 600, 400, 500 77.
400 RCHSZE = RCHSZE + 1 78.
RCHSET(RCHSZE) = NABOR 79.
DEG1 = DEG1 + QSIZE(NABOR) 80.
MARKER(NABOR) = 1 81.
GOTO 600 82.
500 IF ( MARK .GT. 1 ) GOTO 600 83.
NOVRLP = NOVRLP + 1 84.
OVRLP(NOVRLP) = NABOR 85.
MARKER(NABOR) = 2 86.
600 CONTINUE 87.

从交集 中,根据条件( ),确定不可区分的节点。

HEAD 是不可区分节点其 QLINK 中最后一个节点应该设的值,初始为 0 ,MRGSZE 为可以合并的所有不可区分节点数目,初始为 0 。

DO 1100 循环遍历 NOVRLP 中的每一个节点,设为 NODE 。通过 DO 800 循环遍历其所有相邻节点,确定其是否符合条件 ,如果符合条件则其所有相邻节点的MARKER 值必须均不为0 (C1 中的节点有的为1(本次需要消元的节点),有的为 -1(其它);R1 的 MARKER 值为 1,这两部分是由 QMDRCH 设置的。C2 的 MARKER 值为 -1,R2 中一部分为 1 ( R2 - R1 ),一部分为 2 ( ),这两部分是本脚本设置的。除了这些节点以外,其它节点的 MARKER 值并未改变,也就是说,其它节点如何合并进别的节点中了则 MARKER 值为 -1 ,否则为 0)。

如果存在 MARKER 值等于0的相邻节点(这个节点不在 中),则将 NODE 的 MARKER 值从 2 修改为 1(重置为 R1 的默认值), 说明 NODE 不符合条件,跳出 DO 800 循环,GOTO 1100,遍历下一个接待你 。

如果均不存在 MARKER 值等于0的相邻节点,说明 NODE 符合条件,接着向下走,MRGSZE 加上 NODE 的 QSIZE 值。将 NODE 节点的 MARKER 值改为 -1 。

后面的几句脚本 (110-114句) 通过一个 UNTIL 循环找到 NODE 在 QLINK 中的最后一个节点 LNODE。最后 QLINK(LNODE) = HEAD 将 LNODE 的 QLINK 值设为 HEAD ,然后将 HEAD 设为 NODE 。这段脚本的含义是将所有不可区分的节点合并到一块,对于第一个不可区分的节点 NODE,将其QLINK 链接中最后一个节点的 QLINK 值设为 0 (HEAD 初始值);对于第二个不可区分的节点,将其QLINK 链接中最后一个节点的 QLINK 值设为第一个节点,从而将第二个不可区分的节点和第一个不可区分的节点链接起来,以此类推。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
C           --------------------------------------------                  88.
C FROM THE OVERLAPPED SET, DETERMINE THE NODES 89.
C THAT CAN BE MERGED TOGETHER. 90.
C -------------------------------------------- 91.
700 HEAD = 0 92.
MRGSZE = 0 93.
DO 1100 IOV = 1, NOVRLP 94.
NODE = OVRLP(IOV) 95.
JSTRT = XADJ(NODE) 96.
JSTOP = XADJ(NODE+1) - 1 97.
DO 800 J = JSTRT, JSTOP 98.
NABOR = ADJNCY(J) 99.
IF ( MARKER(NABOR) .NE. 0 ) GOTO 800 100.
MARKER(NODE) = 1 101.
GOTO 1100 102.
800 CONTINUE 103.
C ----------------------------------------- 104.
C NODE BELONGS TO THE NEW MERGED SUPERNODE. 105.
C UPDATE THE VECTORS QLINK AND QSIZE. 106.
C ----------------------------------------- 107.
MRGSZE = MRGSZE + QSIZE(NODE) 108.
MARKER(NODE) = - 1 109.
LNODE = NODE 110.
900 LINK = QLINK(LNODE) 111.
IF ( LINK .LE. 0 ) GOTO 1000 112.
LNODE = LINK 113.
GOTO 900 114.
1000 QLINK(LNODE) = HEAD 115.
HEAD = NODE 116.
1100 CONTINUE 117.

如果 HEAD 等于0,那说明符合条件 的节点数目是 0 个,因此没有需要合并的不可区分的节点,那么直接 GOTO 1200

1
IF ( HEAD .LE. 0 )  GOTO 1200                                118.

如果 HEAD 不等于0,那么就说明存在需要合并的不可区分的节点,形成一个 supernode ,这个 supernode 的代表节点就是 HEAD。修改 HEAD 的 QSIZE 为 MRGSZE(所有不可区分节点数目)。

将 HEAD 的 DEG 值设为 DEG0 + DEG1 - 1 ,这里的度数是 。这句话解释如下,这些节点在 中 ,其相邻节点满足条件 ,因此其尚未消元的相邻节点均在 中。进一步消除 之后, 中的节点彼此相连,Head 节点也在 中,因此其与 中其它节点均相连,数目为 DEG0 - 1(减去 HEAD 自己)。而在 中,除了与 共同的节点, 额外与 HEAD 相邻的节点均在 R2 - R1 中,其节点数目为 DEG1 。因此,合并这两点,其 DEG 值更新为 DEG0 + DEG1 - 1

这些节点的 MARKER 值已经设为 -1 了,最后 supernode 的代表节点 HEAD 的 MARKER 值更新为 2 。

1
2
3
QSIZE(HEAD) = MRGSZE                                      119.
DEG(HEAD) = DEG0 + DEG1 - 1 120.
MARKER(HEAD) = 2 121.

最后重置 MARKER 值,将 NBRHD 的节点 ROOT 的 MARKER 值重新设为 0 ,将 RCHSET 中节点 (R2-R1) 的 MARKER 值均重新设为 0 。

1
2
3
4
5
6
7
8
9
10
11
12
13
C           --------------------                                         122.
C RESET MARKER VALUES. 123.
C -------------------- 124.
1200 ROOT = NBRHD(INHD) 125.
MARKER(ROOT) = 0 126.
IF ( RCHSZE .LE. 0 ) GOTO 1400 127.
DO 1300 IRCH = 1, RCHSZE 128.
NODE = RCHSET(IRCH) 129.
MARKER(NODE) = 0 130.
1300 CONTINUE 131.
1400 CONTINUE 132.
RETURN 133.
END 134.

最后除了 NOVRLP 中节点的 MARKER 值修改了之外(不可区分的supernode 为 2,被合并的节点为 -1),其它节点的 MARKER 值均已经还原重置。

注意 NBRHD 中所以所有节点的 MARKER 值均为 0

参考文献

  1. George A, Liu J, Ng E. Computer solution of sparse linear systems[J]. Oak Ridge National Laboratory, 1994.
  • 版权声明: 本博客所有文章除特别声明外,著作权归作者所有。转载请注明出处!
  • Copyrights © 2019-2022 Vincere Zhou
  • 访问人数: | 浏览次数:

请我喝杯茶吧~

支付宝
微信