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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1007 - (hide annotations)
Mon Mar 5 01:00:53 2007 UTC (15 years, 6 months ago) by gross
File MIME type: text/plain
File size: 7564 byte(s)
addSystemMatrix compiles
1 jgs 150 /*
2 elspeth 616 ************************************************************
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 jgs 150 */
12 jgs 82
13     /**************************************************************/
14    
15     /* Finley: Converting an element list into a matrix shape */
16    
17     /**************************************************************/
18    
19 jgs 150 /* Author: gross@access.edu.au */
20     /* Version: $Id$ */
21 jgs 82
22     /**************************************************************/
23    
24     #include "IndexList.h"
25    
26 ksteube 989 /* Translate from distributed/local array indices to global indices */
27    
28     #ifdef PASO_MPI
29 gross 1007 int Finley_IndexList_localToGlobal(Finley_NodeDistribution *dofDistribution, int localIndex) {
30 ksteube 989 /*
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 gross 1007 index_t my_CPU=dofDistribution->MPIInfo->rank;
37 ksteube 989 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 jgs 82 /**************************************************************/
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 ksteube 989 void Finley_IndexList_insertElements(Finley_IndexList* index_list, Finley_Mesh* mesh, Finley_ElementFile* elements,
53 jgs 123 bool_t reduce_row_order, index_t* row_Label,
54     bool_t reduce_col_order, index_t* col_Label) {
55 ksteube 989 /* 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 jgs 123 dim_t e,kr,kc,NN_row,NN_col,i,icol,irow;
58 ksteube 989 #ifdef PASO_MPI
59     num_CPUs = mesh->MPIInfo->size;
60     my_CPU = mesh->MPIInfo->rank;
61     #endif
62     /* print_mesh_statistics( mesh, TRUE ); */
63 jgs 82
64     if (elements!=NULL) {
65 jgs 123 dim_t NN=elements->ReferenceElement->Type->numNodes;
66     index_t id[NN],*row_node,*col_node;
67 jgs 82 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 ksteube 989 if (num_CPUs == 1) {
83 jgs 123 for (color=elements->minColor;color<=elements->maxColor;color++) {
84 jgs 102 #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 jgs 82 }
93 jgs 102 }
94 jgs 82 }
95 jgs 102 }
96     }
97 ksteube 989 }
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 gross 1007 icol = Finley_IndexList_localToGlobal(col_degreeOfFreedomDistribution, icol);
125 ksteube 989 Finley_IndexList_insertIndex(&(index_list[irow]),icol);
126     }
127     }
128     }
129     }
130     }
131     #endif
132     } /* More than one CPU */
133 jgs 82 }
134     return;
135     }
136    
137     /* inserts row index row into the Finley_IndexList in if it does not exist */
138    
139 jgs 123 void Finley_IndexList_insertIndex(Finley_IndexList* in, index_t index) {
140     dim_t i;
141 jgs 82 /* is index in in? */
142     for (i=0;i<in->n;i++) {
143 jgs 102 if (in->index[i]==index) return;
144 jgs 82 }
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 jgs 102 in->extension=TMPMEMALLOC(1,Finley_IndexList);
150 jgs 82 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 jgs 123 dim_t Finley_IndexList_count(Finley_IndexList* in) {
165 jgs 82 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 jgs 123 void Finley_IndexList_toArray(Finley_IndexList* in, index_t* array) {
175     dim_t i;
176 jgs 82 if (in!=NULL) {
177 gross 416 for (i=0;i<in->n;i++) array[i]=in->index[i];
178 jgs 82 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 jgs 150 * 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 jgs 123 * 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 jgs 102 * Revision 1.4 2004/12/15 07:08:32 jgs
203 jgs 97 * *** empty log message ***
204 jgs 123 * Revision 1.1.1.1.2.3 2005/06/29 02:34:50 gross
205     * some changes towards 64 integers in finley
206 jgs 82 *
207 jgs 123 * 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 jgs 97 *
210 jgs 82 *
211 jgs 123 *
212 jgs 82 */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26