/[escript]/trunk/dudley/src/IndexList.cpp
ViewVC logotype

Contents of /trunk/dudley/src/IndexList.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3911 - (show annotations)
Thu Jun 14 01:01:03 2012 UTC (7 years, 3 months ago) by jfenwick
Original Path: trunk/dudley/src/IndexList.c
File MIME type: text/plain
File size: 7300 byte(s)
Copyright changes
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2012 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14 /**************************************************************/
15
16 /* Dudley: Converting an element list into a matrix shape */
17
18 /**************************************************************/
19
20 #include "IndexList.h"
21
22 /* Translate from distributed/local array indices to global indices */
23
24 /**************************************************************/
25 /* inserts the contributions from the element matrices of elements
26 into the row index col. If symmetric is set, only the upper
27 triangle of the matrix is stored. */
28
29 void Dudley_IndexList_insertElements(Dudley_IndexList * index_list, Dudley_ElementFile * elements,
30 bool_t reduce_row_order, index_t * row_map,
31 bool_t reduce_col_order, index_t * col_map)
32 {
33 /* index_list is an array of linked lists. Each entry is a row (DOF) and contains the indices to the non-zero columns */
34 index_t color;
35 dim_t e, kr, kc, NN_row, NN_col, icol, irow, NN;
36 if (elements != NULL)
37 {
38 NN = elements->numNodes;
39 NN_col = (elements->numShapes);
40 NN_row = (elements->numShapes);
41
42 for (color = elements->minColor; color <= elements->maxColor; color++)
43 {
44 #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
45 for (e = 0; e < elements->numElements; e++)
46 {
47 if (elements->Color[e] == color)
48 {
49 for (kr = 0; kr < NN_row; kr++)
50 {
51 irow = row_map[elements->Nodes[INDEX2(kr, e, NN)]];
52 for (kc = 0; kc < NN_col; kc++)
53 {
54 icol = col_map[elements->Nodes[INDEX2(kc, e, NN)]];
55 Dudley_IndexList_insertIndex(&(index_list[irow]), icol);
56 }
57 }
58 }
59 }
60 }
61 }
62 return;
63 }
64
65 void Dudley_IndexList_insertElementsWithRowRange(Dudley_IndexList * index_list, index_t firstRow, index_t lastRow,
66 Dudley_ElementFile * elements, index_t * row_map, index_t * col_map)
67 {
68 /* this does not resolve macro elements */
69 index_t color;
70 dim_t e, kr, kc, icol, irow, NN;
71 if (elements != NULL)
72 {
73 NN = elements->numNodes;
74 for (color = elements->minColor; color <= elements->maxColor; color++)
75 {
76 #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
77 for (e = 0; e < elements->numElements; e++)
78 {
79 if (elements->Color[e] == color)
80 {
81 for (kr = 0; kr < NN; kr++)
82 {
83 irow = row_map[elements->Nodes[INDEX2(kr, e, NN)]];
84 if ((firstRow <= irow) && (irow < lastRow))
85 {
86 irow -= firstRow;
87 for (kc = 0; kc < NN; kc++)
88 {
89 icol = col_map[elements->Nodes[INDEX2(kc, e, NN)]];
90 Dudley_IndexList_insertIndex(&(index_list[irow]), icol);
91 }
92 }
93 }
94 }
95 }
96 }
97 }
98 }
99
100 void Dudley_IndexList_insertElementsWithRowRangeNoMainDiagonal(Dudley_IndexList * index_list, index_t firstRow,
101 index_t lastRow, Dudley_ElementFile * elements,
102 index_t * row_map, index_t * col_map)
103 {
104 /* this does not resolve macro elements */
105 index_t color;
106 dim_t e, kr, kc, icol, irow, NN, irow_loc;
107 if (elements != NULL)
108 {
109 NN = elements->numNodes;
110 for (color = elements->minColor; color <= elements->maxColor; color++)
111 {
112 #pragma omp for private(e,irow,kr,kc,icol,irow_loc) schedule(static)
113 for (e = 0; e < elements->numElements; e++)
114 {
115 if (elements->Color[e] == color)
116 {
117 for (kr = 0; kr < NN; kr++)
118 {
119 irow = row_map[elements->Nodes[INDEX2(kr, e, NN)]];
120 if ((firstRow <= irow) && (irow < lastRow))
121 {
122 irow_loc = irow - firstRow;
123 for (kc = 0; kc < NN; kc++)
124 {
125 icol = col_map[elements->Nodes[INDEX2(kc, e, NN)]];
126 if (icol != irow)
127 Dudley_IndexList_insertIndex(&(index_list[irow_loc]), icol);
128 }
129 }
130 }
131 }
132 }
133 }
134 }
135 }
136
137 /* inserts row index row into the Dudley_IndexList in if it does not exist */
138
139 void Dudley_IndexList_insertIndex(Dudley_IndexList * in, index_t index)
140 {
141 dim_t i;
142 /* is index in in? */
143 for (i = 0; i < in->n; i++)
144 {
145 if (in->index[i] == index)
146 return;
147 }
148 /* index could not be found */
149 if (in->n == INDEXLIST_LENGTH)
150 {
151 /* if in->index is full check the extension */
152 if (in->extension == NULL)
153 {
154 in->extension = TMPMEMALLOC(1, Dudley_IndexList);
155 if (Dudley_checkPtr(in->extension))
156 return;
157 in->extension->n = 0;
158 in->extension->extension = NULL;
159 }
160 Dudley_IndexList_insertIndex(in->extension, index);
161 }
162 else
163 {
164 /* insert index into in->index */
165 in->index[in->n] = index;
166 in->n++;
167 }
168 }
169
170 /* counts the number of row indices in the Dudley_IndexList in */
171
172 dim_t Dudley_IndexList_count(Dudley_IndexList * in, index_t range_min, index_t range_max)
173 {
174 dim_t i;
175 dim_t out = 0;
176 register index_t itmp;
177 if (in == NULL)
178 {
179 return 0;
180 }
181 else
182 {
183 for (i = 0; i < in->n; i++)
184 {
185 itmp = in->index[i];
186 if ((itmp >= range_min) && (range_max > itmp))
187 ++out;
188 }
189 return out + Dudley_IndexList_count(in->extension, range_min, range_max);
190 }
191 }
192
193 /* count the number of row indices in the Dudley_IndexList in */
194
195 void Dudley_IndexList_toArray(Dudley_IndexList * in, index_t * array, index_t range_min, index_t range_max,
196 index_t index_offset)
197 {
198 dim_t i, ptr;
199 register index_t itmp;
200 if (in != NULL)
201 {
202 ptr = 0;
203 for (i = 0; i < in->n; i++)
204 {
205 itmp = in->index[i];
206 if ((itmp >= range_min) && (range_max > itmp))
207 {
208 array[ptr] = itmp + index_offset;
209 ptr++;
210 }
211
212 }
213 Dudley_IndexList_toArray(in->extension, &(array[ptr]), range_min, range_max, index_offset);
214 }
215 }
216
217 /* deallocates the Dudley_IndexList in by recursive calls */
218
219 void Dudley_IndexList_free(Dudley_IndexList * in)
220 {
221 if (in != NULL)
222 {
223 Dudley_IndexList_free(in->extension);
224 TMPMEMFREE(in);
225 }
226 }
227
228 /* creates a Paso_pattern from a range of indices */
229 Paso_Pattern *Dudley_IndexList_createPattern(dim_t n0, dim_t n, Dudley_IndexList * index_list, index_t range_min,
230 index_t range_max, index_t index_offset)
231 {
232 dim_t *ptr = NULL;
233 register dim_t s, i, itmp;
234 index_t *index = NULL;
235 Paso_Pattern *out = NULL;
236
237 ptr = MEMALLOC(n + 1 - n0, index_t);
238 if (!Dudley_checkPtr(ptr))
239 {
240 /* get the number of connections per row */
241 #pragma omp parallel for schedule(static) private(i)
242 for (i = n0; i < n; ++i)
243 {
244 ptr[i - n0] = Dudley_IndexList_count(&index_list[i], range_min, range_max);
245 }
246 /* accumulate ptr */
247 s = 0;
248 for (i = n0; i < n; ++i)
249 {
250 itmp = ptr[i - n0];
251 ptr[i - n0] = s;
252 s += itmp;
253 }
254 ptr[n - n0] = s;
255 /* fill index */
256 index = MEMALLOC(ptr[n - n0], index_t);
257 if (!Dudley_checkPtr(index))
258 {
259 #pragma omp parallel for schedule(static)
260 for (i = n0; i < n; ++i)
261 {
262 Dudley_IndexList_toArray(&index_list[i], &index[ptr[i - n0]], range_min, range_max, index_offset);
263 }
264 out = Paso_Pattern_alloc(MATRIX_FORMAT_DEFAULT, n - n0, range_max + index_offset, ptr, index);
265 }
266 }
267 if (!Dudley_noError())
268 {
269 MEMFREE(ptr);
270 MEMFREE(index);
271 Paso_Pattern_free(out);
272 }
273 return out;
274 }

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26