| 1 |
|
| 2 |
/* $Id$ */ |
| 3 |
|
| 4 |
/******************************************************* |
| 5 |
* |
| 6 |
* Copyright 2003-2007 by ACceSS MNRF |
| 7 |
* Copyright 2007 by University of Queensland |
| 8 |
* |
| 9 |
* http://esscc.uq.edu.au |
| 10 |
* Primary Business: Queensland, Australia |
| 11 |
* Licensed under the Open Software License version 3.0 |
| 12 |
* http://www.opensource.org/licenses/osl-3.0.php |
| 13 |
* |
| 14 |
*******************************************************/ |
| 15 |
|
| 16 |
/**************************************************************/ |
| 17 |
|
| 18 |
/* Finley: Converting an element list into a matrix shape */ |
| 19 |
|
| 20 |
/**************************************************************/ |
| 21 |
|
| 22 |
#include "IndexList.h" |
| 23 |
|
| 24 |
/* Translate from distributed/local array indices to global indices */ |
| 25 |
|
| 26 |
/**************************************************************/ |
| 27 |
/* inserts the contributions from the element matrices of elements |
| 28 |
into the row index col. If symmetric is set, only the upper |
| 29 |
triangle of the matrix is stored. */ |
| 30 |
|
| 31 |
void Finley_IndexList_insertElements(Finley_IndexList* index_list, Finley_ElementFile* elements, |
| 32 |
bool_t reduce_row_order, index_t* row_map, |
| 33 |
bool_t reduce_col_order, index_t* col_map) { |
| 34 |
/* index_list is an array of linked lists. Each entry is a row (DOF) and contains the indices to the non-zero columns */ |
| 35 |
index_t color, *id=NULL; |
| 36 |
dim_t e,kr,kc,NN_row,NN_col,i,icol,irow, NN, *row_node=NULL,*col_node=NULL; |
| 37 |
if (elements!=NULL) { |
| 38 |
NN=elements->numNodes; |
| 39 |
id=TMPMEMALLOC(NN, index_t); |
| 40 |
if (! Finley_checkPtr(id) ) { |
| 41 |
for (i=0;i<NN;i++) id[i]=i; |
| 42 |
if (reduce_col_order) { |
| 43 |
col_node=elements->ReferenceElement->Type->linearNodes; |
| 44 |
NN_col=elements->LinearReferenceElement->Type->numNodes; |
| 45 |
} else { |
| 46 |
col_node=id; |
| 47 |
NN_col=elements->ReferenceElement->Type->numNodes; |
| 48 |
} |
| 49 |
if (reduce_row_order) { |
| 50 |
row_node=elements->ReferenceElement->Type->linearNodes; |
| 51 |
NN_row=elements->LinearReferenceElement->Type->numNodes; |
| 52 |
} else { |
| 53 |
row_node=id; |
| 54 |
NN_row=elements->ReferenceElement->Type->numNodes; |
| 55 |
} |
| 56 |
for (color=elements->minColor;color<=elements->maxColor;color++) { |
| 57 |
#pragma omp for private(e,irow,kr,kc,icol) schedule(static) |
| 58 |
for (e=0;e<elements->numElements;e++) { |
| 59 |
if (elements->Color[e]==color) { |
| 60 |
for (kr=0;kr<NN_row;kr++) { |
| 61 |
irow=row_map[elements->Nodes[INDEX2(row_node[kr],e,NN)]]; |
| 62 |
for (kc=0;kc<NN_col;kc++) { |
| 63 |
icol=col_map[elements->Nodes[INDEX2(col_node[kc],e,NN)]]; |
| 64 |
Finley_IndexList_insertIndex(&(index_list[irow]),icol); |
| 65 |
} |
| 66 |
} |
| 67 |
} |
| 68 |
} |
| 69 |
} |
| 70 |
TMPMEMFREE(id); |
| 71 |
} |
| 72 |
} |
| 73 |
return; |
| 74 |
} |
| 75 |
|
| 76 |
void Finley_IndexList_insertElementsWithRowRange(Finley_IndexList* index_list, index_t firstRow, index_t lastRow, |
| 77 |
Finley_ElementFile* elements, index_t* row_map, index_t* col_map) |
| 78 |
{ |
| 79 |
index_t color; |
| 80 |
dim_t e,kr,kc,icol,irow, NN; |
| 81 |
if (elements!=NULL) { |
| 82 |
NN=elements->numNodes; |
| 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;kr++) { |
| 88 |
irow=row_map[elements->Nodes[INDEX2(kr,e,NN)]]; |
| 89 |
if ((firstRow<=irow) && (irow < lastRow)) { |
| 90 |
irow-=firstRow; |
| 91 |
for (kc=0;kc<NN;kc++) { |
| 92 |
icol=col_map[elements->Nodes[INDEX2(kc,e,NN)]]; |
| 93 |
Finley_IndexList_insertIndex(&(index_list[irow]),icol); |
| 94 |
} |
| 95 |
} |
| 96 |
} |
| 97 |
} |
| 98 |
} |
| 99 |
} |
| 100 |
} |
| 101 |
} |
| 102 |
|
| 103 |
/* inserts row index row into the Finley_IndexList in if it does not exist */ |
| 104 |
|
| 105 |
void Finley_IndexList_insertIndex(Finley_IndexList* in, index_t index) { |
| 106 |
dim_t i; |
| 107 |
/* is index in in? */ |
| 108 |
for (i=0;i<in->n;i++) { |
| 109 |
if (in->index[i]==index) return; |
| 110 |
} |
| 111 |
/* index could not be found */ |
| 112 |
if (in->n==INDEXLIST_LENGTH) { |
| 113 |
/* if in->index is full check the extension */ |
| 114 |
if (in->extension==NULL) { |
| 115 |
in->extension=TMPMEMALLOC(1,Finley_IndexList); |
| 116 |
if (Finley_checkPtr(in->extension)) return; |
| 117 |
in->extension->n=0; |
| 118 |
in->extension->extension=NULL; |
| 119 |
} |
| 120 |
Finley_IndexList_insertIndex(in->extension,index); |
| 121 |
} else { |
| 122 |
/* insert index into in->index*/ |
| 123 |
in->index[in->n]=index; |
| 124 |
in->n++; |
| 125 |
} |
| 126 |
} |
| 127 |
|
| 128 |
/* counts the number of row indices in the Finley_IndexList in */ |
| 129 |
|
| 130 |
dim_t Finley_IndexList_count(Finley_IndexList* in, index_t range_min,index_t range_max) { |
| 131 |
dim_t i; |
| 132 |
dim_t out=0; |
| 133 |
register index_t itmp; |
| 134 |
if (in==NULL) { |
| 135 |
return 0; |
| 136 |
} else { |
| 137 |
for (i=0;i<in->n;i++) { |
| 138 |
itmp=in->index[i]; |
| 139 |
if ((itmp>=range_min) && (range_max>itmp)) ++out; |
| 140 |
} |
| 141 |
return out+Finley_IndexList_count(in->extension, range_min,range_max); |
| 142 |
} |
| 143 |
} |
| 144 |
|
| 145 |
/* count the number of row indices in the Finley_IndexList in */ |
| 146 |
|
| 147 |
void Finley_IndexList_toArray(Finley_IndexList* in, index_t* array, index_t range_min,index_t range_max, index_t index_offset) { |
| 148 |
dim_t i, ptr; |
| 149 |
register index_t itmp; |
| 150 |
if (in!=NULL) { |
| 151 |
ptr=0; |
| 152 |
for (i=0;i<in->n;i++) { |
| 153 |
itmp=in->index[i]; |
| 154 |
if ((itmp>=range_min) && (range_max>itmp)) { |
| 155 |
array[ptr]=itmp+index_offset; |
| 156 |
ptr++; |
| 157 |
} |
| 158 |
|
| 159 |
} |
| 160 |
Finley_IndexList_toArray(in->extension,&(array[ptr]), range_min, range_max, index_offset); |
| 161 |
} |
| 162 |
} |
| 163 |
|
| 164 |
/* deallocates the Finley_IndexList in by recursive calls */ |
| 165 |
|
| 166 |
void Finley_IndexList_free(Finley_IndexList* in) { |
| 167 |
if (in!=NULL) { |
| 168 |
Finley_IndexList_free(in->extension); |
| 169 |
TMPMEMFREE(in); |
| 170 |
} |
| 171 |
} |
| 172 |
|
| 173 |
/* creates a Paso_pattern from a range of indices */ |
| 174 |
Paso_Pattern* Finley_IndexList_createPattern(dim_t n,Finley_IndexList* index_list,index_t range_min,index_t range_max,index_t index_offset) |
| 175 |
{ |
| 176 |
dim_t *ptr=NULL; |
| 177 |
register dim_t s,i,itmp; |
| 178 |
index_t *index=NULL; |
| 179 |
Paso_Pattern* out=NULL; |
| 180 |
|
| 181 |
ptr=MEMALLOC(n+1,index_t); |
| 182 |
if (! Finley_checkPtr(ptr) ) { |
| 183 |
/* get the number of connections per row */ |
| 184 |
#pragma omp parallel for schedule(static) private(i) |
| 185 |
for(i=0;i<n;++i) { |
| 186 |
ptr[i]=Finley_IndexList_count(&index_list[i],range_min,range_max); |
| 187 |
} |
| 188 |
/* accumulate ptr */ |
| 189 |
s=0; |
| 190 |
for(i=0;i<n;++i) { |
| 191 |
itmp=ptr[i]; |
| 192 |
ptr[i]=s; |
| 193 |
s+=itmp; |
| 194 |
} |
| 195 |
ptr[n]=s; |
| 196 |
/* fill index */ |
| 197 |
index=MEMALLOC(ptr[n],index_t); |
| 198 |
if (! Finley_checkPtr(index)) { |
| 199 |
#pragma omp parallel for schedule(static) |
| 200 |
for(i=0;i<n;++i) { |
| 201 |
Finley_IndexList_toArray(&index_list[i],&index[ptr[i]],range_min,range_max,index_offset); |
| 202 |
} |
| 203 |
out=Paso_Pattern_alloc(PATTERN_FORMAT_DEFAULT,1,1,n,ptr,index); |
| 204 |
} |
| 205 |
} |
| 206 |
if (! Finley_noError()) { |
| 207 |
MEMFREE(ptr); |
| 208 |
MEMFREE(index); |
| 209 |
Paso_Pattern_free(out); |
| 210 |
} |
| 211 |
return out; |
| 212 |
} |
Parent Directory
| 
Revision Log