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

Contents of /trunk/finley/src/IndexList.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 969 - (show annotations)
Tue Feb 13 23:02:23 2007 UTC (12 years, 9 months ago) by ksteube
File MIME type: text/plain
File size: 7640 byte(s)
Parallelization using MPI for solution of implicit problems.

Parallelization for explicit problems has already been accomplished in
the main SVN branch.

This is incomplete and is not ready for use.


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 my_CPU, 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 if (localIndex < dofDistribution->numLocal) {
37 localIndex = localIndex + dofDistribution->vtxdist[my_CPU];
38 }
39 else {
40 localIndex = dofDistribution->indexExternal[localIndex-dofDistribution->numLocal];
41 }
42 return(localIndex);
43 }
44 #endif
45
46 /**************************************************************/
47 /* inserts the contributions from the element matrices of elements
48 into the row index col. If symmetric is set, only the upper
49 triangle of the matrix is stored. */
50
51 void Finley_IndexList_insertElements(Finley_IndexList* index_list, Finley_Mesh* mesh, Finley_ElementFile* elements,
52 bool_t reduce_row_order, index_t* row_Label,
53 bool_t reduce_col_order, index_t* col_Label) {
54 /* index_list is an array of linked lists. Each entry is a row (DOF) and contains the indices to the non-zero columns */
55 index_t color, num_CPUs = 1, my_CPU = 0;
56 dim_t e,kr,kc,NN_row,NN_col,i,icol,irow;
57 #ifdef PASO_MPI
58 num_CPUs = mesh->MPIInfo->size;
59 my_CPU = mesh->MPIInfo->rank;
60 #endif
61 /* print_mesh_statistics( mesh, TRUE ); */
62
63 if (elements!=NULL) {
64 dim_t NN=elements->ReferenceElement->Type->numNodes;
65 index_t id[NN],*row_node,*col_node;
66 for (i=0;i<NN;i++) id[i]=i;
67 if (reduce_col_order) {
68 col_node=elements->ReferenceElement->Type->linearNodes;
69 NN_col=elements->LinearReferenceElement->Type->numNodes;
70 } else {
71 col_node=id;
72 NN_col=elements->ReferenceElement->Type->numNodes;
73 }
74 if (reduce_row_order) {
75 row_node=elements->ReferenceElement->Type->linearNodes;
76 NN_row=elements->LinearReferenceElement->Type->numNodes;
77 } else {
78 row_node=id;
79 NN_row=elements->ReferenceElement->Type->numNodes;
80 }
81 if (num_CPUs == 1) {
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_row;kr++) {
87 irow=row_Label[elements->Nodes[INDEX2(row_node[kr],e,NN)]];
88 for (kc=0;kc<NN_col;kc++) {
89 icol=col_Label[elements->Nodes[INDEX2(col_node[kc],e,NN)]];
90 Finley_IndexList_insertIndex(&(index_list[irow]),icol);
91 }
92 }
93 }
94 }
95 }
96 }
97 else { /* More than one CPU (what's below should also work for one CPU, but let's build confidence in it first) */
98 #ifdef PASO_MPI
99 Finley_NodeDistribution *row_degreeOfFreedomDistribution;
100 Finley_NodeDistribution *col_degreeOfFreedomDistribution;
101 if (reduce_col_order) {
102 col_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
103 }
104 else {
105 col_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
106 }
107 if (reduce_row_order) {
108 row_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
109 }
110 else {
111 row_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
112 }
113 /* Not using loop over colors as above */ {
114 #pragma omp for private(e,irow,kr,kc,icol) schedule(static)
115 for (e=0;e<elements->numElements;e++) {
116 for (kr=0;kr<NN_row;kr++) {
117 irow=row_Label[elements->Nodes[INDEX2(row_node[kr],e,NN)]];
118 if (irow < row_degreeOfFreedomDistribution->numLocal) {
119 for (kc=0;kc<NN_col;kc++) {
120 /* Get the local col ID */
121 icol=col_Label[elements->Nodes[INDEX2(col_node[kc],e,NN)]];
122 /* Convert to global col ID (row ID is saved as local value) */
123 icol = Finley_IndexList_localToGlobal(col_degreeOfFreedomDistribution, my_CPU, icol);
124 Finley_IndexList_insertIndex(&(index_list[irow]),icol);
125 printf("ksteube Finley_IndexList_insertIndex cpu= %d irow= %d icol= %d\n", my_CPU, 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