/[escript]/trunk-mpi-branch/finley/src/IndexList.c
ViewVC logotype

Contents of /trunk-mpi-branch/finley/src/IndexList.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1007 - (show annotations)
Mon Mar 5 01:00:53 2007 UTC (13 years, 8 months ago) by gross
File MIME type: text/plain
File size: 7564 byte(s)
addSystemMatrix compiles
1 /*
2 ************************************************************
3 * Copyright 2006 by ACcESS MNRF *
4 * *
5 * http://www.access.edu.au *
6 * Primary Business: Queensland, Australia *
7 * Licensed under the Open Software License version 3.0 *
8 * http://www.opensource.org/licenses/osl-3.0.php *
9 * *
10 ************************************************************
11 */
12
13 /**************************************************************/
14
15 /* Finley: Converting an element list into a matrix shape */
16
17 /**************************************************************/
18
19 /* Author: gross@access.edu.au */
20 /* Version: $Id$ */
21
22 /**************************************************************/
23
24 #include "IndexList.h"
25
26 /* Translate from distributed/local array indices to global indices */
27
28 #ifdef PASO_MPI
29 int Finley_IndexList_localToGlobal(Finley_NodeDistribution *dofDistribution, int localIndex) {
30 /*
31 get global id of icol
32 if icol is internal node (on this CPU): use icol+vtxdist[my_CPU]
33 else use indexExternal[icol-numLocal] to get global index of node
34 (actually DOF...the NodeDistribution structure should have been called DofDistribution)
35 */
36 index_t my_CPU=dofDistribution->MPIInfo->rank;
37 if (localIndex < dofDistribution->numLocal) {
38 localIndex = localIndex + dofDistribution->vtxdist[my_CPU];
39 }
40 else {
41 localIndex = dofDistribution->indexExternal[localIndex-dofDistribution->numLocal];
42 }
43 return(localIndex);
44 }
45 #endif
46
47 /**************************************************************/
48 /* inserts the contributions from the element matrices of elements
49 into the row index col. If symmetric is set, only the upper
50 triangle of the matrix is stored. */
51
52 void Finley_IndexList_insertElements(Finley_IndexList* index_list, Finley_Mesh* mesh, Finley_ElementFile* elements,
53 bool_t reduce_row_order, index_t* row_Label,
54 bool_t reduce_col_order, index_t* col_Label) {
55 /* index_list is an array of linked lists. Each entry is a row (DOF) and contains the indices to the non-zero columns */
56 index_t color, num_CPUs = 1, my_CPU = 0;
57 dim_t e,kr,kc,NN_row,NN_col,i,icol,irow;
58 #ifdef PASO_MPI
59 num_CPUs = mesh->MPIInfo->size;
60 my_CPU = mesh->MPIInfo->rank;
61 #endif
62 /* print_mesh_statistics( mesh, TRUE ); */
63
64 if (elements!=NULL) {
65 dim_t NN=elements->ReferenceElement->Type->numNodes;
66 index_t id[NN],*row_node,*col_node;
67 for (i=0;i<NN;i++) id[i]=i;
68 if (reduce_col_order) {
69 col_node=elements->ReferenceElement->Type->linearNodes;
70 NN_col=elements->LinearReferenceElement->Type->numNodes;
71 } else {
72 col_node=id;
73 NN_col=elements->ReferenceElement->Type->numNodes;
74 }
75 if (reduce_row_order) {
76 row_node=elements->ReferenceElement->Type->linearNodes;
77 NN_row=elements->LinearReferenceElement->Type->numNodes;
78 } else {
79 row_node=id;
80 NN_row=elements->ReferenceElement->Type->numNodes;
81 }
82 if (num_CPUs == 1) {
83 for (color=elements->minColor;color<=elements->maxColor;color++) {
84 #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
85 for (e=0;e<elements->numElements;e++) {
86 if (elements->Color[e]==color) {
87 for (kr=0;kr<NN_row;kr++) {
88 irow=row_Label[elements->Nodes[INDEX2(row_node[kr],e,NN)]];
89 for (kc=0;kc<NN_col;kc++) {
90 icol=col_Label[elements->Nodes[INDEX2(col_node[kc],e,NN)]];
91 Finley_IndexList_insertIndex(&(index_list[irow]),icol);
92 }
93 }
94 }
95 }
96 }
97 }
98 else { /* More than one CPU (what's below should also work for one CPU, but let's build confidence in it first) */
99 #ifdef PASO_MPI
100 Finley_NodeDistribution *row_degreeOfFreedomDistribution;
101 Finley_NodeDistribution *col_degreeOfFreedomDistribution;
102 if (reduce_col_order) {
103 col_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
104 }
105 else {
106 col_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
107 }
108 if (reduce_row_order) {
109 row_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
110 }
111 else {
112 row_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
113 }
114 /* Not using loop over colors as above */ {
115 #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
116 for (e=0;e<elements->numElements;e++) {
117 for (kr=0;kr<NN_row;kr++) {
118 irow=row_Label[elements->Nodes[INDEX2(row_node[kr],e,NN)]];
119 if (irow < row_degreeOfFreedomDistribution->numLocal) {
120 for (kc=0;kc<NN_col;kc++) {
121 /* Get the local col ID */
122 icol=col_Label[elements->Nodes[INDEX2(col_node[kc],e,NN)]];
123 /* Convert to global col ID (row ID is saved as local value) */
124 icol = Finley_IndexList_localToGlobal(col_degreeOfFreedomDistribution, icol);
125 Finley_IndexList_insertIndex(&(index_list[irow]),icol);
126 }
127 }
128 }
129 }
130 }
131 #endif
132 } /* More than one CPU */
133 }
134 return;
135 }
136
137 /* inserts row index row into the Finley_IndexList in if it does not exist */
138
139 void Finley_IndexList_insertIndex(Finley_IndexList* in, index_t index) {
140 dim_t i;
141 /* is index in in? */
142 for (i=0;i<in->n;i++) {
143 if (in->index[i]==index) return;
144 }
145 /* index could not be found */
146 if (in->n==INDEXLIST_LENGTH) {
147 /* if in->index is full check the extension */
148 if (in->extension==NULL) {
149 in->extension=TMPMEMALLOC(1,Finley_IndexList);
150 if (Finley_checkPtr(in->extension)) return;
151 in->extension->n=0;
152 in->extension->extension=NULL;
153 }
154 Finley_IndexList_insertIndex(in->extension,index);
155 } else {
156 /* insert index into in->index*/
157 in->index[in->n]=index;
158 in->n++;
159 }
160 }
161
162 /* counts the number of row indices in the Finley_IndexList in */
163
164 dim_t Finley_IndexList_count(Finley_IndexList* in) {
165 if (in==NULL) {
166 return 0;
167 } else {
168 return (in->n)+Finley_IndexList_count(in->extension);
169 }
170 }
171
172 /* count the number of row indices in the Finley_IndexList in */
173
174 void Finley_IndexList_toArray(Finley_IndexList* in, index_t* array) {
175 dim_t i;
176 if (in!=NULL) {
177 for (i=0;i<in->n;i++) array[i]=in->index[i];
178 Finley_IndexList_toArray(in->extension,&(array[in->n]));
179 }
180 }
181
182 /* deallocates the Finley_IndexList in by recursive calls */
183
184 void Finley_IndexList_free(Finley_IndexList* in) {
185 if (in!=NULL) {
186 Finley_IndexList_free(in->extension);
187 TMPMEMFREE(in);
188 }
189 }
190
191 /*
192 * $Log$
193 * Revision 1.6 2005/09/15 03:44:22 jgs
194 * Merge of development branch dev-02 back to main trunk on 2005-09-15
195 *
196 * Revision 1.5.2.1 2005/09/07 06:26:18 gross
197 * the solver from finley are put into the standalone package paso now
198 *
199 * Revision 1.5 2005/07/08 04:07:51 jgs
200 * Merge of development branch back to main trunk on 2005-07-08
201 *
202 * Revision 1.4 2004/12/15 07:08:32 jgs
203 * *** empty log message ***
204 * Revision 1.1.1.1.2.3 2005/06/29 02:34:50 gross
205 * some changes towards 64 integers in finley
206 *
207 * Revision 1.1.1.1.2.2 2004/11/24 01:37:13 gross
208 * some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now
209 *
210 *
211 *
212 */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26