# Annotation of /trunk/paso/src/Pattern_reduceBandwidth.cpp

Revision 4803 - (hide annotations)
Wed Mar 26 06:52:28 2014 UTC (5 years, 11 months ago) by caltinay
File size: 10210 byte(s)
```Removed obsolete wrappers for malloc and friends.
Paso_Pattern -> paso::Pattern

```
 1 ksteube 1315 2 jfenwick 3981 /***************************************************************************** 3 ksteube 1811 * 4 jfenwick 4657 * Copyright (c) 2003-2014 by University of Queensland 5 jfenwick 3981 6 ksteube 1811 * 7 * Primary Business: Queensland, Australia 8 * Licensed under the Open Software License version 3.0 9 10 * 11 jfenwick 3981 * Development until 2012 by Earth Systems Science Computational Center (ESSCC) 12 jfenwick 4657 * Development 2012-2013 by School of Earth Sciences 13 * Development from 2014 by Centre for Geoscience Computing (GeoComp) 14 jfenwick 3981 * 15 *****************************************************************************/ 16 ksteube 1315 17 ksteube 1811 18 jfenwick 3981 /************************************************************************************/ 19 ksteube 1315 20 caltinay 3642 /* Paso: pattern: a simple algorithm to create a new labeling */ 21 /* of input/output with a minimum bandwidth. The algorithm */ 22 /* starts from a vertex with a minimum number of neighbours */ 23 /* (connections in the pattern). From this root it searches */ 24 /* recursively the neighbours. In the last level a new root */ 25 /* is picked (again the vertex with the minimum number of */ 26 /* neighbours) and the process is repeated. The algorithm is */ 27 /* repeated until the maximum number of vertices in a level */ 28 /* cannot be reduced. */ 29 ksteube 1315 30 jfenwick 3981 /************************************************************************************/ 31 ksteube 1315 32 jfenwick 2608 /* Author: Lutz Gross, l.gross@uq.edu.au */ 33 ksteube 1315 34 jfenwick 3981 /************************************************************************************/ 35 ksteube 1315 36 37 #include "Paso.h" 38 #include "Pattern.h" 39 #include "Common.h" 40 41 caltinay 4803 namespace paso { 42 43 caltinay 3642 /* calculate initial bandwidth for a given labeling */ 44 caltinay 4803 dim_t Pattern_getBandwidth(Pattern* pattern, index_t* label) { 45 ksteube 1315 register index_t k; 46 index_t iptr; 47 dim_t bandwidth = 0, local_bandwidth=0,i ; 48 49 #pragma omp parallel private(local_bandwidth) 50 { 51 local_bandwidth=0; 52 #pragma omp for private(i,iptr,k) 53 for (i=0;inumOutput;++i) { 54 k=label[i]; 55 for (iptr=pattern->ptr[i];iptrptr[i+1];++iptr) { 56 local_bandwidth = MAX(local_bandwidth, ABS(k - label[pattern->index[iptr]])); 57 } 58 } 59 #pragma omp critical 60 { 61 bandwidth=MAX(local_bandwidth,bandwidth); 62 } 63 } 64 return bandwidth; 65 } 66 /* structures used to create an ordering by increasing degree */ 67 typedef struct Paso_DegreeAndIdx { 68 dim_t deg; 69 index_t idx; 70 } Paso_DegreeAndIdx; 71 72 int Paso_comparDegreeAndIdx(const void *arg1,const void *arg2){ 73 Paso_DegreeAndIdx a1=*(Paso_DegreeAndIdx*)arg1; 74 Paso_DegreeAndIdx a2=*(Paso_DegreeAndIdx*)arg2; 75 if (a1.dega2.deg) { 78 return 1; 79 } else if (a1.idxa2.idx) { 82 return 1; 83 } else { 84 return 0; 85 } 86 } 87 88 caltinay 4803 /* Pattern_dropTree drops a tree in pattern from root */ 89 ksteube 1315 90 /* root - on input the starting point of the tree. */ 91 /* AssignedLevel- array of length pattern->numOutput indicating the level assigned to vertex */ 92 /* VerticesInTree - on output contains the vertices used in tree (array of length pattern->numOutput) */ 93 caltinay 3642 /* numLevels- on output the number of levels used. */ 94 ksteube 1315 /* firstVertexInLevel - on output firstVertexInLevel[i] points to the first vertex of level i in VerticesInTree. */ 95 /* firstVertexInLevel[i+1]-firstVertexInLevel[i] is the number of vertices in level i. */ 96 /* (array of length pattern->numOutput+1) */ 97 /* max_LevelWidth_abort- input param which triggers early return if LevelWidth becomes >= max_LevelWidth_abort */ 98 caltinay 4803 bool Pattern_dropTree(index_t root, 99 Pattern *pattern, 100 ksteube 1315 index_t *AssignedLevel, 101 index_t *VerticesInTree, 102 dim_t* numLevels, 103 index_t* firstVertexInLevel, 104 dim_t max_LevelWidth_abort, 105 caltinay 4803 dim_t N) 106 ksteube 1315 { 107 dim_t nlvls,i; 108 index_t level_top,k ,itest,j; 109 110 #pragma omp parallel for private(i) 111 for (i=0;inumInput;++i) AssignedLevel[i]=-1; 112 113 nlvls=0; 114 AssignedLevel[root] = 0; 115 VerticesInTree[0] = root; 116 firstVertexInLevel[0]=0; 117 level_top=firstVertexInLevel[0]+1; 118 while (firstVertexInLevel[nlvls]=max_LevelWidth_abort) return FALSE; 122 for (i=firstVertexInLevel[nlvls-1];iptr[k];jptr[k+1];++j) { 125 itest = pattern->index[j]; 126 if (AssignedLevel[itest]<0) { 127 #ifdef BOUNDS_CHECK 128 caltinay 4803 if (itest < 0 || itest >= N) { printf("BOUNDS_CHECK %s %d itest=%d\n", __FILE__, __LINE__, itest); exit(1); } 129 if (level_top < 0 || level_top >= N) { printf("BOUNDS_CHECK %s %d level_top=%d\n", __FILE__, __LINE__, level_top); exit(1); } 130 ksteube 1315 #endif 131 AssignedLevel[itest] = nlvls; 132 VerticesInTree[level_top] = itest; 133 level_top++; 134 } 135 } 136 } 137 } 138 *numLevels=nlvls; 139 return TRUE; 140 } 141 142 /* the driver */ 143 caltinay 4803 void Pattern_reduceBandwidth(Pattern* pattern,index_t* oldToNew) { 144 jfenwick 1981 dim_t i, initial_bandwidth, N=pattern->numOutput, numLabledVertices, numLevels, max_LevelWidth, min_deg,deg, numVerticesInTree=0, bandwidth; 145 ksteube 1315 Paso_DegreeAndIdx* degAndIdx=NULL; 146 index_t root, *AssignedLevel=NULL, *VerticesInTree=NULL, *firstVertexInLevel=NULL,k, *oldLabel=NULL; 147 /* check input */ 148 if (N != pattern->numInput) { 149 caltinay 4803 Esys_setError(VALUE_ERROR,"Pattern_reduceBandwidth: pattern needs to be for a square matrix."); 150 caltinay 2873 } else if (N > 0) { 151 gross 1841 /* printf("relabeling of %d DOFs started.\n",N); */ 152 jfenwick 4291 degAndIdx=new Paso_DegreeAndIdx[N]; 153 jfenwick 4521 oldLabel=new index_t[N]; 154 jfenwick 4291 AssignedLevel=new index_t[N]; 155 VerticesInTree=new index_t[N]; 156 firstVertexInLevel=new index_t[N+1]; 157 jfenwick 3259 if (! ( Esys_checkPtr(degAndIdx) || Esys_checkPtr(oldLabel) || Esys_checkPtr(AssignedLevel) || Esys_checkPtr(VerticesInTree) || Esys_checkPtr(firstVertexInLevel) ) ) { 158 ksteube 1315 /* get the initial bandwidth */ 159 #pragma omp parallel for private(i) 160 for (i=0;iptr[i+1]-pattern->ptr[i]; 169 } 170 /* create an ordering with increasing degree */ 171 #ifdef USE_QSORTG 172 qsortG(degAndIdx,(size_t)N,sizeof(Paso_DegreeAndIdx),Paso_comparDegreeAndIdx); 173 #else 174 qsort(degAndIdx,(size_t)N,sizeof(Paso_DegreeAndIdx),Paso_comparDegreeAndIdx); 175 #endif 176 177 root=degAndIdx[0].idx; 178 numLabledVertices=0; 179 180 while (root>=0) { 181 max_LevelWidth=N+1; 182 183 caltinay 4803 while (Pattern_dropTree(root,pattern,AssignedLevel,VerticesInTree, 184 ksteube 1315 &numLevels,firstVertexInLevel,max_LevelWidth,N) ) { 185 186 /* find new maximum level width */ 187 max_LevelWidth=0; 188 for (i=0;i= N+1) { printf("BOUNDS_CHECK %s %d i=%d N=%d\n", __FILE__, __LINE__, i, N); exit(1); } 191 ksteube 1315 #endif 192 max_LevelWidth=MAX(max_LevelWidth,firstVertexInLevel[i+1]-firstVertexInLevel[i]); 193 } 194 /* find a vertex in the last level which has minimum degree */ 195 min_deg=N+1; 196 root=-1; 197 for (i = firstVertexInLevel[numLevels-1]; iptr[k+1]-pattern->ptr[k]; 200 if (deg= N) { printf("BOUNDS_CHECK %s %d i=%d numLabeledVertices=%d root=%d N=%d firstVertexInLevel[numLevels]=%d\n", __FILE__, __LINE__, i, numLabledVertices, root, N, firstVertexInLevel[numLevels]); exit(1); } 210 ksteube 1315 #endif 211 oldLabel[numLabledVertices+i]=VerticesInTree[i]; 212 } 213 } 214 /* now the vertices in the current tree */ 215 for (i=0;i= N) { printf("BOUNDS_CHECK %s %d i=%d numLabeledVertices=%d root=%d N=%d\n", __FILE__, __LINE__, i, numLabledVertices, root, N); exit(1); } 218 ksteube 1315 #endif 219 oldToNew[oldLabel[numLabledVertices+i]]=numLabledVertices+i; 220 } 221 numLabledVertices+=numVerticesInTree; 222 caltinay 3642 /* new search for a vertex which is not labeled yet */ 223 ksteube 1315 root=-1; 224 for (i=0;i=initial_bandwidth) { 234 gross 1841 /* printf("initial labeling used.\n"); */ 235 ksteube 1315 #pragma omp parallel for private(i) 236 for (i=0;i

## Properties

Name Value
svn:mergeinfo /branches/amg_from_3530/paso/src/Pattern_reduceBandwidth.cpp:3531-3826 /branches/lapack2681/paso/src/Pattern_reduceBandwidth.cpp:2682-2741 /branches/pasowrap/paso/src/Pattern_reduceBandwidth.cpp:3661-3674 /branches/py3_attempt2/paso/src/Pattern_reduceBandwidth.cpp:3871-3891 /branches/restext/paso/src/Pattern_reduceBandwidth.cpp:2610-2624 /branches/ripleygmg_from_3668/paso/src/Pattern_reduceBandwidth.cpp:3669-3791 /branches/stage3.0/paso/src/Pattern_reduceBandwidth.cpp:2569-2590 /branches/symbolic_from_3470/paso/src/Pattern_reduceBandwidth.cpp:3471-3974 /branches/symbolic_from_3470/ripley/test/python/paso/src/Pattern_reduceBandwidth.cpp:3517-3974 /release/3.0/paso/src/Pattern_reduceBandwidth.cpp:2591-2601 /trunk/paso/src/Pattern_reduceBandwidth.cpp:4257-4344 /trunk/ripley/test/python/paso/src/Pattern_reduceBandwidth.cpp:3480-3515