| Diploma Thesis Percolation Simulation C++ Sourcecode Documentation |
00001 // cellcount.h 00002 // the core ClusterCounter algorithm 00003 // counts best if only applied to small cells of spheres 00004 // 00005 // (c) 2000 Andreas Krueger 00006 // Version 1.0 00007 // last change: 28.6.2000 00008 00009 // the naive cell counters... 00010 00011 namespace counters{ 00012 00013 // this is the inner loop of the NOW-used counter 00014 void neighbour_recursion (sphere *array, NUMBER first, NUMBER last, NUMBER actualsph, NUMBER clnumber) { 00015 00016 array[actualsph].clno=clnumber; // set this sphere as "counted" by assigning clusternumber 00017 00018 sphere temp( getdim(array[actualsph].c) ); // used as a temporary distance vector 00019 00020 for (NUMBER nextsph=first;nextsph<=last;nextsph++) { 00021 if ( array[nextsph].clno==0 ) { 00022 if ( overlap2(array[actualsph],array[nextsph], temp) ) { 00023 // this overlap-function uses a temp-sphere-coordinatevector 00024 // for storing the difference of the two vectors 00025 // (it is much faster than creating one temporary vector 00026 // for each overlap-test by using overlap (sphere, sphere) ) 00027 neighbour_recursion (array, first, last, nextsph, clnumber); 00028 } 00029 } 00030 } 00031 } 00032 00033 // used NOW 00034 NUMBER count_clusters_recursively (sphere *array, 00035 NUMBER first, 00036 NUMBER last, 00037 NUMBER first_clno) 00038 { 00039 NUMBER clusternumber=first_clno; 00040 for (NUMBER k=first; k<=last; k++) { 00041 if (array[k].clno==0) { // if yet uncounted 00042 neighbour_recursion(array, k,last, k, clusternumber); // the actual counting process 00043 clusternumber++; // increase cluster-counter 00044 } 00045 } 00046 return (clusternumber); 00047 } 00048 00049 // this naive recursion is used NOW (per divide&conquer in find_ffc, into_file8 ) 00050 NUMBER list_to_temp_array_and_find_cluster(sphere *spheres, 00051 NUMLIST &sphlist, 00052 NUMBER first_clno) 00053 { 00054 00055 sphere *copy= new sphere[sphlist.size()+1]; 00056 00057 NUMLIST::iterator iter; 00058 NUMBER index=0; 00059 for (iter = sphlist.begin(); iter != sphlist.end(); iter++){ 00060 index++; 00061 copy[index]=spheres[*iter]; 00062 } 00063 00064 NUMBER first_free_clusternumber=count_clusters_recursively(copy,1,index,first_clno); 00065 00066 index=0; 00067 for (iter = sphlist.begin(); iter != sphlist.end(); iter++){ 00068 index++; 00069 spheres[*iter]=copy[index]; 00070 } 00071 00072 delete[] copy; 00073 00074 return (first_free_clusternumber); 00075 } 00076 00077 00078 NUMBER neighbour_naive (sphere *array, NUMBER first, NUMBER last, NUMBER actualsph, NUMBER clnumber) { 00079 array[actualsph].clno=clnumber; 00080 00081 myVector distance(getdim(array[actualsph].c)); 00082 bool overlapped; 00083 00084 NUMBER size=1; 00085 00086 for (NUMBER nextsph=first;nextsph<=last;nextsph++) { 00087 if ( array[nextsph].clno==0 ) { 00088 00089 // if ( overlap(array[actualsph],array[nextsph]) ) { 00090 // the following 3 lines are faster: 00091 // because of not using operator-(vect1, vect2) which needs a vector temp; 00092 00093 distance=(array[nextsph].c); 00094 distance-=(array[actualsph].c); // using -= instead of - is much faster (no vector temp) 00095 overlapped=length(distance) < (array[actualsph].r + array[nextsph].r) ; 00096 00097 if (overlapped) { 00098 00099 size+=neighbour_naive (array, first, last, nextsph, clnumber); // RECURSION! 00100 } 00101 } 00102 } 00103 return size; 00104 } 00105 00106 clock_t countclusters_naive(sphere *array, NUMBER first, NUMBER last, NUMBER* clustertable) { 00107 clock_t start=clock(); 00108 00109 NUMBER size; 00110 NUMBER clusternumber=0; // start with # 1 00111 00112 for (NUMBER k=first; k<=last; k++) { 00113 if (array[k].clno==0) { // if yet uncounted 00114 00115 clusternumber++; // increase cluster-counter 00116 size=neighbour_naive(array, k,last, k, clusternumber); // the actual counting process 00117 clustertable[clusternumber]=size; // store size to table of all clusters 00118 00119 } 00120 } 00121 clustertable[clusternumber+1]=0; // mark the end in that table by: "last cluster is of size 0" 00122 return (clock()-start); 00123 } 00124 00125 00126 // the naive cell counters for lists of spheres instead of arrays... 00127 // (unfortunately slower) 00128 00129 NUMBER neighbour_l (sphere *array, 00130 NUMLIST &sphlist, // search only in this part of the array 00131 NUMLIST::iterator firstsph, // from firstsph to sphlist.end() 00132 NUMLIST::iterator actualsph, // search for overlap with this sphere 00133 NUMBER clnumber) 00134 { 00135 array[*actualsph].clno=clnumber; // this sphere is set to be found in this cluster 00136 00137 myVector distance(getdim(array[*actualsph].c)); 00138 bool overlapped; 00139 00140 NUMBER size=1; 00141 00142 NUMLIST::iterator nextsph; 00143 00144 for (nextsph = firstsph; nextsph != sphlist.end(); nextsph++){ 00145 00146 if ( array[*nextsph].clno==0 ) { 00147 00148 //if ( overlap2(array[*actualsph],array[*nextsph]) ) { 00149 // the following 4 lines do the same, but are faster than that: 00150 // because of not using operator-(vect1, vect2) which needs a myVector temp; 00151 00152 distance=(array[*nextsph].c); 00153 distance-=(array[*actualsph].c); 00154 overlapped=length(distance) < (array[*actualsph].r + array[*nextsph].r) ; 00155 if (overlapped) { 00156 size+=neighbour_l (array, sphlist, firstsph, nextsph, clnumber); // recursion, now with nextsph 00157 } 00158 } 00159 } 00160 return size; 00161 } 00162 00163 00164 clock_t countclusters_l(sphere *array, 00165 NUMLIST &sphlist, // only search the part of the array specified by sphlist 00166 NUMBER *clustertable, // returns this part of the clustertable (ended by 0) 00167 NUMBER &clusternumber) // and the next clusternumber 00168 00169 { 00170 clock_t start=clock(); 00171 00172 NUMBER size; 00173 00174 NUMLIST::iterator k; 00175 for (k = sphlist.begin(); k != sphlist.end(); k++){ // search for neighbours from each sphere 00176 00177 if (array[*k].clno==0) { // if yet uncounted 00178 00179 size=neighbour_l(array,sphlist, k, k, clusternumber); // the actual counting process 00180 clustertable[clusternumber]=size; // store size to table of all clusters 00181 clusternumber++; // increase cluster-counter 00182 } 00183 } 00184 clustertable[clusternumber]=0; // mark the end in that table by: "last cluster is of size 0" 00185 return (clock()-start); // return the time (and the clustertable and the next clusternumber!) 00186 } 00187 00188 00189 NUMBER list_to_given_array_and_find_cluster(sphere *spheres, 00190 sphere *copy, 00191 NUMLIST &sphlist, 00192 NUMBER first_clno) 00193 { 00194 NUMLIST::iterator iter; 00195 NUMBER index=0; 00196 for (iter = sphlist.begin(); iter != sphlist.end(); iter++){ 00197 index++; 00198 copy[index]=spheres[*iter]; 00199 } 00200 00201 NUMBER first_free_clusternumber=count_clusters_recursively(copy,1,index,first_clno); 00202 00203 index=0; 00204 for (iter = sphlist.begin(); iter != sphlist.end(); iter++){ 00205 index++; 00206 spheres[*iter]=copy[index]; 00207 } 00208 return (first_free_clusternumber); 00209 } 00210 00211 00212 NUMBER count_clusters_iteratively (sphere *array, 00213 NUMBER first, 00214 NUMBER last, 00215 NUMBER first_clno) 00216 { 00217 // the idea to this algorithm can be found in 00218 // S.D.Stoddard, J.Comp.Phys 27, 291-293 (1978) 00219 // Identifying clusters in computer experiments on systems of particles 00220 00221 NUMBER i,j,k; 00222 NUMBER N=last-first+1; 00223 NUMBER* neighbour=new NUMBER[N]; // afterwards contains a "list" of linked spheres 00224 for (i=0;i<N;i++) {neighbour[i]=i;} // initialize link"list" 00225 NUMBER ntemp; // used for swapping 00226 sphere temp( getdim(array[first].c) ); // used as a temporary distance vector 00227 00228 for (i=0; i<N-1; i++) { 00229 if (neighbour[i]==i) { 00230 j=i; 00231 do { 00232 for (k=i+1; k<N; k++) { 00233 if (neighbour[k]==k){ 00234 if (overlap2(array[first+j],array[first+k], temp)){ 00235 ntemp=neighbour[j]; 00236 neighbour[j]=neighbour[k]; 00237 neighbour[k]=ntemp; 00238 } 00239 } 00240 } 00241 j=neighbour[j]; 00242 } while (j!=i); 00243 } 00244 } 00245 // cout<<"\n"; 00246 // for (i=0;i<N;i++){cout<<"[L("<<i<<")="<<neighbour[i]<<"] ";} 00247 00248 NUMBER clusternumber=first_clno; 00249 for (i=0; i<N; i++){ 00250 if (array[first+i].clno==0){ 00251 j=i; 00252 do { 00253 array[first+i].clno=clusternumber; 00254 i=neighbour[i]; 00255 } while (i!=j); 00256 clusternumber++; 00257 } 00258 } 00259 00260 delete[] neighbour; 00261 return clusternumber; 00262 } 00263 00264 NUMBER count_clusters_iteratively2 (sphere *array, 00265 NUMBER first, 00266 NUMBER last, 00267 NUMBER first_clno) 00268 { 00269 // the idea to this algorithm can be found in 00270 // S.D.Stoddard, J.Comp.Phys 27, 291-293 (1978) 00271 // Identifying clusters in computer experiments on systems of particles 00272 00273 NUMBER i,j,k; 00274 NUMBER N=last-first+1; 00275 00276 NUMBER* neighbour=new NUMBER[N]; // first~ 1 ; last~ last-first+1 00277 for (i=0;i<N;i++) {neighbour[i]=i;} // initialize linklist 00278 00279 NUMBER ntemp; // used for swapping 00280 sphere temp( getdim(array[first].c) ); // used as a temporary distance vector 00281 00282 NUMBER clusternumber=first_clno; 00283 00284 for (i=0; i<N-1; i++) { 00285 if (neighbour[i]==i) { 00286 array[first+i].clno=clusternumber; 00287 j=i; 00288 do { 00289 for (k=i+1; k<N; k++) { 00290 if (neighbour[k]==k){ 00291 if (overlap2(array[first+j],array[first+k], temp)){ 00292 array[first+k].clno=clusternumber; 00293 ntemp=neighbour[j]; 00294 neighbour[j]=neighbour[k]; 00295 neighbour[k]=ntemp; 00296 } 00297 } 00298 } 00299 j=neighbour[j]; 00300 } while (j!=i); 00301 clusternumber++; 00302 } 00303 } 00304 00305 delete [] neighbour; 00306 return clusternumber; 00307 } 00308 00309 NUMBER count_clusters_iteratively_only_selected_spherenumbers (sphere *array, 00310 std::vector<NUMBER> &sphnumber, 00311 NUMBER N, 00312 NUMBER first_clno) 00313 { 00314 // the idea to this algorithm can be found in 00315 // S.D.Stoddard, J.Comp.Phys 27, 291-293 (1978) 00316 // Identifying clusters in computer experiments on systems of particles 00317 00318 NUMBER i,j,k; 00319 NUMBER* neighbour=new NUMBER[N+1]; // afterwards contains a "list" of linked spheres 00320 for (i=1;i<=N;i++) {neighbour[i]=i;} // initialize link"list" 00321 NUMBER ntemp; // used for swapping 00322 sphere temp( getdim(array[sphnumber[N]].c) ); // used as a temporary distance vector 00323 00324 for (i=1; i<=N-1; i++) { 00325 if (neighbour[i]==i) { 00326 j=i; 00327 do { 00328 for (k=i+1; k<=N; k++) { 00329 if (neighbour[k]==k){ 00330 if (overlap2(array[sphnumber[j]],array[sphnumber[k]], temp)){ 00331 ntemp=neighbour[j]; 00332 neighbour[j]=neighbour[k]; 00333 neighbour[k]=ntemp; 00334 } 00335 } 00336 } 00337 j=neighbour[j]; 00338 } while (j!=i); 00339 } 00340 } 00341 // cout<<"\n"; 00342 // for (i=0;i<N;i++){cout<<"[L("<<i<<")="<<neighbour[i]<<"] ";} 00343 00344 NUMBER clusternumber=first_clno; 00345 for (i=1; i<=N; i++){ 00346 if (array[sphnumber[i]].clno==0){ 00347 j=i; 00348 do { 00349 array[sphnumber[i]].clno=clusternumber; 00350 i=neighbour[i]; 00351 } while (i!=j); 00352 clusternumber++; 00353 } 00354 } 00355 00356 delete[] neighbour; 00357 return clusternumber; 00358 } 00359 00360 00361 00362 } // end of namespace counters
| Diploma Thesis Sourcecode
Documentation check out the text and the executable binaries |