Newer
Older
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
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
SUBROUTINE GRAPH_SLOAN(N,NE,INPN,NPN,XNPN,IADJ,ADJ,XADJ)
************************************************************************
*
* PURPOSE:
* --------
*
* Form adjacency list for a graph corresponding to a finite element
* mesh
*
* INPUT:
* -----
*
* N - Number of nodes in graph (finite element mesh)
* NE - Number of elements in finite element mesh
* INPN - Length of NPN = XNPN(NE+1)-1
* NPN - List of node numbers for each element
* XNPN - Index vector for NPN
* - nodes for element I are found in NPN(J), where
* J = XNPN(I), XNPN(I)+1, ..., XNPN(I+1)-1
* IADJ - Length of vector ADJ
* - Set IADJ=NE*NEN*(NEN-1) for a mesh of a single type of
* element with NEN nodes
* - IADJ=(NEN(1)*(NEN(1)-1)+,.....,+NEN(NE)*(NEN(NE)-1))
* for a mesh of elements with varying numbers of nodes
* ADJ - Undefined
* XADJ - Undefined
*
* OUTPUT:
* -------
*
* N - Unchanged
* NE - Unchanged
* INPN - Unchanged
* NPN - Unchanged
* XNPN - Unchanged
* IADJ - Unchanged
* ADJ - Adjacency list for all nodes in graph
* - List of length 2E where E is the number of edges in
* the graph (note that 2E = XADJ(N+1)-1 )
* XADJ - Index vector for ADJ
* - Nodes adjacent to node I are found in ADJ(J), where
* J = XADJ(I), XADJ(I)+1, ..., XADJ(I+1)-1
* - Degree of node I given by XADJ(I+1)-XADJ(I)
*
* NOTES:
* ------
*
* This routine typically requires about 25 percent elbow room for
* assembling the ADJ list (i.e. IADJ/2E is typically around 1.25).
* In some cases, the elbow room may be larger (IADJ/2E is slightly
* less than 2 for the 3-noded triangle) and in other cases it may be
* zero (IADJ/2E = 1 for bar elements)
*
* PROGRAMMER: Scott Sloan
* -----------
*
* LAST MODIFIED: 1 March 1991 Scott Sloan
* --------------
*
* COPYRIGHT 1989: Scott Sloan
* --------------- Department of Civil Engineering
* University of Newcastle
* NSW 2308
*
************************************************************************
INTEGER I,J,K,L,M,N
INTEGER NE
INTEGER IADJ,INPN,NEN1
INTEGER JSTOP,JSTRT,LSTOP,LSTRT,MSTOP,MSTRT,NODEJ,NODEK
INTEGER ADJ(IADJ),NPN(INPN)
INTEGER XADJ(N+1),XNPN(NE+1)
*
* Initialise the adjacency list and its index vector
*
DO 5 I=1,IADJ
ADJ(I)=0
5 CONTINUE
DO 10 I=1,N
XADJ(I)=0
10 CONTINUE
*
* Estimate the degree of each node (always an overestimate)
*
DO 30 I=1,NE
JSTRT=XNPN(I)
JSTOP=XNPN(I+1)-1
NEN1 =JSTOP-JSTRT
DO 20 J=JSTRT,JSTOP
NODEJ=NPN(J)
XADJ(NODEJ)=XADJ(NODEJ)+NEN1
20 CONTINUE
30 CONTINUE
*
* Reconstruct XADJ to point to start of each set of neighbours
*
L=1
DO 40 I=1,N
L=L+XADJ(I)
XADJ(I)=L-XADJ(I)
40 CONTINUE
XADJ(N+1)=L
*
* Form adjacency list (which may contain zeros)
*
DO 90 I=1,NE
JSTRT=XNPN(I)
JSTOP=XNPN(I+1)-1
DO 80 J=JSTRT,JSTOP-1
NODEJ=NPN(J)
LSTRT=XADJ(NODEJ)
LSTOP=XADJ(NODEJ+1)-1
DO 70 K=J+1,JSTOP
NODEK=NPN(K)
DO 50 L=LSTRT,LSTOP
IF(ADJ(L).EQ.NODEK)GOTO 70
IF(ADJ(L).EQ.0)GOTO 55
50 CONTINUE
WRITE(6,1000)
STOP
55 CONTINUE
ADJ(L)=NODEK
MSTRT=XADJ(NODEK)
MSTOP=XADJ(NODEK+1)-1
DO 60 M=MSTRT,MSTOP
IF(ADJ(M).EQ.0)GOTO 65
60 CONTINUE
WRITE(6,1000)
STOP
65 CONTINUE
ADJ(M)=NODEJ
70 CONTINUE
80 CONTINUE
90 CONTINUE
*
* Strip any zeros from adjacency list
*
K=0
JSTRT=1
DO 110 I=1,N
JSTOP=XADJ(I+1)-1
DO 100 J=JSTRT,JSTOP
IF(ADJ(J).EQ.0)GOTO 105
K=K+1
ADJ(K)=ADJ(J)
100 CONTINUE
105 CONTINUE
XADJ(I+1)=K+1
JSTRT=JSTOP+1
110 CONTINUE
*
* Error message
*
1000 FORMAT(//,1X,'***ERROR IN GRAPH***',
+ //,1X,'CANNOT ASSEMBLE NODE ADJACENCY LIST',
+ //,1X,'CHECK NPN AND XNPN ARRAYS')
END