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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2548 - (hide annotations)
Mon Jul 20 06:20:06 2009 UTC (10 years, 2 months ago) by jfenwick
Original Path: trunk/finley/src/IndexList.c
File MIME type: text/plain
File size: 8458 byte(s)
Updating copyright notices
1 jgs 82
2 ksteube 1312 /*******************************************************
3 ksteube 1811 *
4 jfenwick 2548 * Copyright (c) 2003-2009 by University of Queensland
5 ksteube 1811 * 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 jgs 82
14 ksteube 1811
15 jgs 82 /**************************************************************/
16    
17 ksteube 1312 /* Finley: Converting an element list into a matrix shape */
18 jgs 82
19     /**************************************************************/
20    
21     #include "IndexList.h"
22    
23 ksteube 1312 /* Translate from distributed/local array indices to global indices */
24    
25 jgs 82 /**************************************************************/
26     /* inserts the contributions from the element matrices of elements
27     into the row index col. If symmetric is set, only the upper
28     triangle of the matrix is stored. */
29    
30 ksteube 971 void Finley_IndexList_insertElements(Finley_IndexList* index_list, Finley_ElementFile* elements,
31 ksteube 1312 bool_t reduce_row_order, index_t* row_map,
32     bool_t reduce_col_order, index_t* col_map) {
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 gross 1028 index_t color, *id=NULL;
35     dim_t e,kr,kc,NN_row,NN_col,i,icol,irow, NN, *row_node=NULL,*col_node=NULL;
36 jgs 82 if (elements!=NULL) {
37 ksteube 1312 NN=elements->numNodes;
38 gross 1028 id=TMPMEMALLOC(NN, index_t);
39     if (! Finley_checkPtr(id) ) {
40 ksteube 1312 for (i=0;i<NN;i++) id[i]=i;
41     if (reduce_col_order) {
42     col_node=elements->ReferenceElement->Type->linearNodes;
43     NN_col=elements->LinearReferenceElement->Type->numNodes;
44     } else {
45     col_node=id;
46     NN_col=elements->ReferenceElement->Type->numNodes;
47     }
48     if (reduce_row_order) {
49     row_node=elements->ReferenceElement->Type->linearNodes;
50     NN_row=elements->LinearReferenceElement->Type->numNodes;
51     } else {
52     row_node=id;
53     NN_row=elements->ReferenceElement->Type->numNodes;
54     }
55     for (color=elements->minColor;color<=elements->maxColor;color++) {
56     #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
57     for (e=0;e<elements->numElements;e++) {
58     if (elements->Color[e]==color) {
59     for (kr=0;kr<NN_row;kr++) {
60     irow=row_map[elements->Nodes[INDEX2(row_node[kr],e,NN)]];
61     for (kc=0;kc<NN_col;kc++) {
62     icol=col_map[elements->Nodes[INDEX2(col_node[kc],e,NN)]];
63     Finley_IndexList_insertIndex(&(index_list[irow]),icol);
64     }
65     }
66     }
67     }
68     }
69     TMPMEMFREE(id);
70     }
71 jgs 82 }
72     return;
73     }
74    
75 ksteube 1312 void Finley_IndexList_insertElementsWithRowRange(Finley_IndexList* index_list, index_t firstRow, index_t lastRow,
76     Finley_ElementFile* elements, index_t* row_map, index_t* col_map)
77     {
78     index_t color;
79 phornby 1628 dim_t e,kr,kc,icol,irow, NN;
80 ksteube 1312 if (elements!=NULL) {
81     NN=elements->numNodes;
82     for (color=elements->minColor;color<=elements->maxColor;color++) {
83     #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
84     for (e=0;e<elements->numElements;e++) {
85     if (elements->Color[e]==color) {
86     for (kr=0;kr<NN;kr++) {
87     irow=row_map[elements->Nodes[INDEX2(kr,e,NN)]];
88     if ((firstRow<=irow) && (irow < lastRow)) {
89     irow-=firstRow;
90     for (kc=0;kc<NN;kc++) {
91     icol=col_map[elements->Nodes[INDEX2(kc,e,NN)]];
92     Finley_IndexList_insertIndex(&(index_list[irow]),icol);
93     }
94     }
95     }
96     }
97     }
98     }
99     }
100     }
101 gross 1722 void Finley_IndexList_insertElementsWithRowRangeNoMainDiagonal(Finley_IndexList* index_list, index_t firstRow, index_t lastRow,
102     Finley_ElementFile* elements, index_t* row_map, index_t* col_map)
103     {
104     index_t color;
105     dim_t e,kr,kc,icol,irow, NN,irow_loc;
106     if (elements!=NULL) {
107     NN=elements->numNodes;
108     for (color=elements->minColor;color<=elements->maxColor;color++) {
109     #pragma omp for private(e,irow,kr,kc,icol,irow_loc) schedule(static)
110     for (e=0;e<elements->numElements;e++) {
111     if (elements->Color[e]==color) {
112     for (kr=0;kr<NN;kr++) {
113     irow=row_map[elements->Nodes[INDEX2(kr,e,NN)]];
114     if ((firstRow<=irow) && (irow < lastRow)) {
115     irow_loc=irow-firstRow;
116     for (kc=0;kc<NN;kc++) {
117     icol=col_map[elements->Nodes[INDEX2(kc,e,NN)]];
118     if (icol != irow) Finley_IndexList_insertIndex(&(index_list[irow_loc]),icol);
119     }
120     }
121     }
122     }
123     }
124     }
125     }
126     }
127 ksteube 1312
128 jgs 82 /* inserts row index row into the Finley_IndexList in if it does not exist */
129    
130 jgs 123 void Finley_IndexList_insertIndex(Finley_IndexList* in, index_t index) {
131     dim_t i;
132 jgs 82 /* is index in in? */
133     for (i=0;i<in->n;i++) {
134 jgs 102 if (in->index[i]==index) return;
135 jgs 82 }
136     /* index could not be found */
137     if (in->n==INDEXLIST_LENGTH) {
138     /* if in->index is full check the extension */
139     if (in->extension==NULL) {
140 jgs 102 in->extension=TMPMEMALLOC(1,Finley_IndexList);
141 jgs 82 if (Finley_checkPtr(in->extension)) return;
142     in->extension->n=0;
143     in->extension->extension=NULL;
144     }
145     Finley_IndexList_insertIndex(in->extension,index);
146     } else {
147     /* insert index into in->index*/
148     in->index[in->n]=index;
149     in->n++;
150     }
151     }
152    
153     /* counts the number of row indices in the Finley_IndexList in */
154    
155 ksteube 1312 dim_t Finley_IndexList_count(Finley_IndexList* in, index_t range_min,index_t range_max) {
156     dim_t i;
157     dim_t out=0;
158     register index_t itmp;
159 jgs 82 if (in==NULL) {
160     return 0;
161     } else {
162 ksteube 1312 for (i=0;i<in->n;i++) {
163     itmp=in->index[i];
164     if ((itmp>=range_min) && (range_max>itmp)) ++out;
165     }
166     return out+Finley_IndexList_count(in->extension, range_min,range_max);
167 jgs 82 }
168     }
169    
170     /* count the number of row indices in the Finley_IndexList in */
171    
172 ksteube 1312 void Finley_IndexList_toArray(Finley_IndexList* in, index_t* array, index_t range_min,index_t range_max, index_t index_offset) {
173     dim_t i, ptr;
174     register index_t itmp;
175 jgs 82 if (in!=NULL) {
176 ksteube 1312 ptr=0;
177     for (i=0;i<in->n;i++) {
178     itmp=in->index[i];
179     if ((itmp>=range_min) && (range_max>itmp)) {
180     array[ptr]=itmp+index_offset;
181     ptr++;
182     }
183    
184     }
185     Finley_IndexList_toArray(in->extension,&(array[ptr]), range_min, range_max, index_offset);
186 jgs 82 }
187     }
188    
189     /* deallocates the Finley_IndexList in by recursive calls */
190    
191     void Finley_IndexList_free(Finley_IndexList* in) {
192     if (in!=NULL) {
193     Finley_IndexList_free(in->extension);
194     TMPMEMFREE(in);
195     }
196     }
197    
198 ksteube 1312 /* creates a Paso_pattern from a range of indices */
199 gross 1552 Paso_Pattern* Finley_IndexList_createPattern(dim_t n0, dim_t n,Finley_IndexList* index_list,index_t range_min,index_t range_max,index_t index_offset)
200 ksteube 1312 {
201     dim_t *ptr=NULL;
202     register dim_t s,i,itmp;
203     index_t *index=NULL;
204     Paso_Pattern* out=NULL;
205    
206 gross 1552 ptr=MEMALLOC(n+1-n0,index_t);
207 ksteube 1312 if (! Finley_checkPtr(ptr) ) {
208     /* get the number of connections per row */
209     #pragma omp parallel for schedule(static) private(i)
210 gross 1552 for(i=n0;i<n;++i) {
211     ptr[i-n0]=Finley_IndexList_count(&index_list[i],range_min,range_max);
212 ksteube 1312 }
213     /* accumulate ptr */
214     s=0;
215 gross 1552 for(i=n0;i<n;++i) {
216     itmp=ptr[i-n0];
217     ptr[i-n0]=s;
218 ksteube 1312 s+=itmp;
219     }
220 gross 1552 ptr[n-n0]=s;
221 ksteube 1312 /* fill index */
222 gross 1552 index=MEMALLOC(ptr[n-n0],index_t);
223 ksteube 1312 if (! Finley_checkPtr(index)) {
224     #pragma omp parallel for schedule(static)
225 gross 1552 for(i=n0;i<n;++i) {
226     Finley_IndexList_toArray(&index_list[i],&index[ptr[i-n0]],range_min,range_max,index_offset);
227 ksteube 1312 }
228 gross 1736 out=Paso_Pattern_alloc(PATTERN_FORMAT_DEFAULT,1,1,n-n0,range_max+index_offset,ptr,index);
229 ksteube 1312 }
230     }
231     if (! Finley_noError()) {
232     MEMFREE(ptr);
233     MEMFREE(index);
234     Paso_Pattern_free(out);
235     }
236     return out;
237     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26