1 2 /* 3 This file contains routines for sorting integers. Values are sorted in place. 4 */ 5 #include <petscsys.h> /*I "petscsys.h" I*/ 6 7 #define SWAP(a,b,t) {t=a;a=b;b=t;} 8 9 #define MEDIAN3(v,a,b,c) \ 10 (v[a]<v[b] \ 11 ? (v[b]<v[c] \ 12 ? b \ 13 : (v[a]<v[c] ? c : a)) \ 14 : (v[c]<v[b] \ 15 ? b \ 16 : (v[a]<v[c] ? a : c))) 17 18 #define MEDIAN(v,right) MEDIAN3(v,right/4,right/2,right/4*3) 19 20 /* -----------------------------------------------------------------------*/ 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "PetscSortInt_Private" 24 /* 25 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 26 Assumes 0 origin for v, number of elements = right+1 (right is index of 27 right-most member). 28 */ 29 static void PetscSortInt_Private(PetscInt *v,PetscInt right) 30 { 31 PetscInt i,j,pivot,tmp; 32 33 if (right <= 1) { 34 if (right == 1) { 35 if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 36 } 37 return; 38 } 39 i = MEDIAN(v,right); /* Choose a pivot */ 40 SWAP(v[0],v[i],tmp); /* Move it out of the way */ 41 pivot = v[0]; 42 for (i=0,j=right+1;;) { 43 while (++i < j && v[i] <= pivot); /* Scan from the left */ 44 while (v[--j] > pivot); /* Scan from the right */ 45 if (i >= j) break; 46 SWAP(v[i],v[j],tmp); 47 } 48 SWAP(v[0],v[j],tmp); /* Put pivot back in place. */ 49 PetscSortInt_Private(v,j-1); 50 PetscSortInt_Private(v+j+1,right-(j+1)); 51 } 52 53 #undef __FUNCT__ 54 #define __FUNCT__ "PetscSortInt" 55 /*@ 56 PetscSortInt - Sorts an array of integers in place in increasing order. 57 58 Not Collective 59 60 Input Parameters: 61 + n - number of values 62 - i - array of integers 63 64 Level: intermediate 65 66 Concepts: sorting^ints 67 68 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 69 @*/ 70 PetscErrorCode PetscSortInt(PetscInt n,PetscInt i[]) 71 { 72 PetscInt j,k,tmp,ik; 73 74 PetscFunctionBegin; 75 if (n<8) { 76 for (k=0; k<n; k++) { 77 ik = i[k]; 78 for (j=k+1; j<n; j++) { 79 if (ik > i[j]) { 80 SWAP(i[k],i[j],tmp); 81 ik = i[k]; 82 } 83 } 84 } 85 } else { 86 PetscSortInt_Private(i,n-1); 87 } 88 PetscFunctionReturn(0); 89 } 90 91 #undef __FUNCT__ 92 #define __FUNCT__ "PetscSortRemoveDupsInt" 93 /*@ 94 PetscSortRemoveDupsInt - Sorts an array of integers in place in increasing order removes all duplicate entries 95 96 Not Collective 97 98 Input Parameters: 99 + n - number of values 100 - ii - array of integers 101 102 Output Parameter: 103 . n - number of non-redundant values 104 105 Level: intermediate 106 107 Concepts: sorting^ints 108 109 .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt() 110 @*/ 111 PetscErrorCode PetscSortRemoveDupsInt(PetscInt *n,PetscInt ii[]) 112 { 113 PetscErrorCode ierr; 114 PetscInt i,s = 0,N = *n, b = 0; 115 116 PetscFunctionBegin; 117 ierr = PetscSortInt(N,ii);CHKERRQ(ierr); 118 for (i=0; i<N-1; i++) { 119 if (ii[b+s+1] != ii[b]) {ii[b+1] = ii[b+s+1]; b++;} 120 else s++; 121 } 122 *n = N - s; 123 PetscFunctionReturn(0); 124 } 125 126 127 /* -----------------------------------------------------------------------*/ 128 #define SWAP2(a,b,c,d,t) {t=a;a=b;b=t;t=c;c=d;d=t;} 129 130 #undef __FUNCT__ 131 #define __FUNCT__ "PetscSortIntWithArray_Private" 132 /* 133 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 134 Assumes 0 origin for v, number of elements = right+1 (right is index of 135 right-most member). 136 */ 137 static PetscErrorCode PetscSortIntWithArray_Private(PetscInt *v,PetscInt *V,PetscInt right) 138 { 139 PetscErrorCode ierr; 140 PetscInt i,vl,last,tmp; 141 142 PetscFunctionBegin; 143 if (right <= 1) { 144 if (right == 1) { 145 if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp); 146 } 147 PetscFunctionReturn(0); 148 } 149 SWAP2(v[0],v[right/2],V[0],V[right/2],tmp); 150 vl = v[0]; 151 last = 0; 152 for (i=1; i<=right; i++) { 153 if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);} 154 } 155 SWAP2(v[0],v[last],V[0],V[last],tmp); 156 ierr = PetscSortIntWithArray_Private(v,V,last-1);CHKERRQ(ierr); 157 ierr = PetscSortIntWithArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 158 PetscFunctionReturn(0); 159 } 160 161 #undef __FUNCT__ 162 #define __FUNCT__ "PetscSortIntWithArray" 163 /*@ 164 PetscSortIntWithArray - Sorts an array of integers in place in increasing order; 165 changes a second array to match the sorted first array. 166 167 Not Collective 168 169 Input Parameters: 170 + n - number of values 171 . i - array of integers 172 - I - second array of integers 173 174 Level: intermediate 175 176 Concepts: sorting^ints with array 177 178 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt() 179 @*/ 180 PetscErrorCode PetscSortIntWithArray(PetscInt n,PetscInt i[],PetscInt Ii[]) 181 { 182 PetscErrorCode ierr; 183 PetscInt j,k,tmp,ik; 184 185 PetscFunctionBegin; 186 if (n<8) { 187 for (k=0; k<n; k++) { 188 ik = i[k]; 189 for (j=k+1; j<n; j++) { 190 if (ik > i[j]) { 191 SWAP2(i[k],i[j],Ii[k],Ii[j],tmp); 192 ik = i[k]; 193 } 194 } 195 } 196 } else { 197 ierr = PetscSortIntWithArray_Private(i,Ii,n-1);CHKERRQ(ierr); 198 } 199 PetscFunctionReturn(0); 200 } 201 202 203 #define SWAP3(a,b,c,d,e,f,t) {t=a;a=b;b=t;t=c;c=d;d=t;t=e;e=f;f=t;} 204 205 #undef __FUNCT__ 206 #define __FUNCT__ "PetscSortIntWithArrayPair_Private" 207 /* 208 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 209 Assumes 0 origin for v, number of elements = right+1 (right is index of 210 right-most member). 211 */ 212 static PetscErrorCode PetscSortIntWithArrayPair_Private(PetscInt *L,PetscInt *J, PetscInt *K, PetscInt right) 213 { 214 PetscErrorCode ierr; 215 PetscInt i,vl,last,tmp; 216 217 PetscFunctionBegin; 218 if (right <= 1) { 219 if (right == 1) { 220 if (L[0] > L[1]) SWAP3(L[0],L[1],J[0],J[1],K[0],K[1],tmp); 221 } 222 PetscFunctionReturn(0); 223 } 224 SWAP3(L[0],L[right/2],J[0],J[right/2],K[0],K[right/2],tmp); 225 vl = L[0]; 226 last = 0; 227 for (i=1; i<=right; i++) { 228 if (L[i] < vl) {last++; SWAP3(L[last],L[i],J[last],J[i],K[last],K[i],tmp);} 229 } 230 SWAP3(L[0],L[last],J[0],J[last],K[0],K[last],tmp); 231 ierr = PetscSortIntWithArrayPair_Private(L,J,K,last-1);CHKERRQ(ierr); 232 ierr = PetscSortIntWithArrayPair_Private(L+last+1,J+last+1,K+last+1,right-(last+1));CHKERRQ(ierr); 233 PetscFunctionReturn(0); 234 } 235 236 #undef __FUNCT__ 237 #define __FUNCT__ "PetscSortIntWithArrayPair" 238 /*@ 239 PetscSortIntWithArrayPair - Sorts an array of integers in place in increasing order; 240 changes a pair of integer arrays to match the sorted first array. 241 242 Not Collective 243 244 Input Parameters: 245 + n - number of values 246 . I - array of integers 247 . J - second array of integers (first array of the pair) 248 - K - third array of integers (second array of the pair) 249 250 Level: intermediate 251 252 Concepts: sorting^ints with array pair 253 254 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortIntWithArray() 255 @*/ 256 PetscErrorCode PetscSortIntWithArrayPair(PetscInt n,PetscInt *L,PetscInt *J, PetscInt *K) 257 { 258 PetscErrorCode ierr; 259 PetscInt j,k,tmp,ik; 260 261 PetscFunctionBegin; 262 if (n<8) { 263 for (k=0; k<n; k++) { 264 ik = L[k]; 265 for (j=k+1; j<n; j++) { 266 if (ik > L[j]) { 267 SWAP3(L[k],L[j],J[k],J[j],K[k],K[j],tmp); 268 ik = L[k]; 269 } 270 } 271 } 272 } else { 273 ierr = PetscSortIntWithArrayPair_Private(L,J,K,n-1);CHKERRQ(ierr); 274 } 275 PetscFunctionReturn(0); 276 } 277 278 279 #undef __FUNCT__ 280 #define __FUNCT__ "PetscSortMPIIntWithArray_Private" 281 /* 282 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 283 Assumes 0 origin for v, number of elements = right+1 (right is index of 284 right-most member). 285 */ 286 static PetscErrorCode PetscSortMPIIntWithArray_Private(PetscMPIInt *v,PetscMPIInt *V,PetscMPIInt right) 287 { 288 PetscErrorCode ierr; 289 PetscMPIInt i,vl,last,tmp; 290 291 PetscFunctionBegin; 292 if (right <= 1) { 293 if (right == 1) { 294 if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp); 295 } 296 PetscFunctionReturn(0); 297 } 298 SWAP2(v[0],v[right/2],V[0],V[right/2],tmp); 299 vl = v[0]; 300 last = 0; 301 for (i=1; i<=right; i++) { 302 if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);} 303 } 304 SWAP2(v[0],v[last],V[0],V[last],tmp); 305 ierr = PetscSortMPIIntWithArray_Private(v,V,last-1);CHKERRQ(ierr); 306 ierr = PetscSortMPIIntWithArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 307 PetscFunctionReturn(0); 308 } 309 310 #undef __FUNCT__ 311 #define __FUNCT__ "PetscSortMPIIntWithArray" 312 /*@ 313 PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order; 314 changes a second array to match the sorted first array. 315 316 Not Collective 317 318 Input Parameters: 319 + n - number of values 320 . i - array of integers 321 - I - second array of integers 322 323 Level: intermediate 324 325 Concepts: sorting^ints with array 326 327 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt() 328 @*/ 329 PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt i[],PetscMPIInt Ii[]) 330 { 331 PetscErrorCode ierr; 332 PetscMPIInt j,k,tmp,ik; 333 334 PetscFunctionBegin; 335 if (n<8) { 336 for (k=0; k<n; k++) { 337 ik = i[k]; 338 for (j=k+1; j<n; j++) { 339 if (ik > i[j]) { 340 SWAP2(i[k],i[j],Ii[k],Ii[j],tmp); 341 ik = i[k]; 342 } 343 } 344 } 345 } else { 346 ierr = PetscSortMPIIntWithArray_Private(i,Ii,n-1);CHKERRQ(ierr); 347 } 348 PetscFunctionReturn(0); 349 } 350 351 /* -----------------------------------------------------------------------*/ 352 #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;} 353 354 #undef __FUNCT__ 355 #define __FUNCT__ "PetscSortIntWithScalarArray_Private" 356 /* 357 Modified from PetscSortIntWithArray_Private(). 358 */ 359 static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right) 360 { 361 PetscErrorCode ierr; 362 PetscInt i,vl,last,tmp; 363 PetscScalar stmp; 364 365 PetscFunctionBegin; 366 if (right <= 1) { 367 if (right == 1) { 368 if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp); 369 } 370 PetscFunctionReturn(0); 371 } 372 SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp); 373 vl = v[0]; 374 last = 0; 375 for (i=1; i<=right; i++) { 376 if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);} 377 } 378 SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp); 379 ierr = PetscSortIntWithScalarArray_Private(v,V,last-1);CHKERRQ(ierr); 380 ierr = PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 381 PetscFunctionReturn(0); 382 } 383 384 #undef __FUNCT__ 385 #define __FUNCT__ "PetscSortIntWithScalarArray" 386 /*@ 387 PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order; 388 changes a second SCALAR array to match the sorted first INTEGER array. 389 390 Not Collective 391 392 Input Parameters: 393 + n - number of values 394 . i - array of integers 395 - I - second array of scalars 396 397 Level: intermediate 398 399 Concepts: sorting^ints with array 400 401 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 402 @*/ 403 PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[]) 404 { 405 PetscErrorCode ierr; 406 PetscInt j,k,tmp,ik; 407 PetscScalar stmp; 408 409 PetscFunctionBegin; 410 if (n<8) { 411 for (k=0; k<n; k++) { 412 ik = i[k]; 413 for (j=k+1; j<n; j++) { 414 if (ik > i[j]) { 415 SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp); 416 ik = i[k]; 417 } 418 } 419 } 420 } else { 421 ierr = PetscSortIntWithScalarArray_Private(i,Ii,n-1);CHKERRQ(ierr); 422 } 423 PetscFunctionReturn(0); 424 } 425 426 427 428 #undef __FUNCT__ 429 #define __FUNCT__ "PetscMergeIntArrayPair" 430 /*@ 431 PetscMergeIntArrayPair - Merges two SORTED integer arrays along with an additional array of integers. 432 The additional arrays are the same length as sorted arrays and are merged 433 in the order determined by the merging of the sorted pair. 434 435 Not Collective 436 437 Input Parameters: 438 + an - number of values in the first array 439 . aI - first sorted array of integers 440 . aJ - first additional array of integers 441 . bn - number of values in the second array 442 . bI - second array of integers 443 - bJ - second additional of integers 444 445 Output Parameters: 446 + n - number of values in the merged array (== an + bn) 447 . I - merged sorted array 448 - J - merged additional array 449 450 Level: intermediate 451 452 Concepts: merging^arrays 453 454 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 455 @*/ 456 PetscErrorCode PetscMergeIntArrayPair(PetscInt an,const PetscInt *aI, const PetscInt *aJ, PetscInt bn, const PetscInt *bI, const PetscInt *bJ, PetscInt *n, PetscInt **L, PetscInt **J) 457 { 458 PetscErrorCode ierr; 459 PetscInt n_, *L_ = *L, *J_= *J, ak, bk, k; 460 461 n_ = an + bn; 462 *n = n_; 463 if(!L_) { 464 ierr = PetscMalloc(n_*sizeof(PetscInt), L); CHKERRQ(ierr); 465 L_ = *L; 466 } 467 if(!J_){ 468 ierr = PetscMalloc(n_*sizeof(PetscInt), &J_); CHKERRQ(ierr); 469 J_ = *J; 470 } 471 k = ak = bk = 0; 472 while(ak < an && bk < bn) { 473 if(aI[ak] <= bI[bk]) { 474 L_[k] = aI[ak]; 475 J_[k] = aJ[ak]; 476 ++ak; 477 ++k; 478 } 479 else { 480 L_[k] = bI[bk]; 481 J_[k] = bJ[bk]; 482 ++bk; 483 ++k; 484 } 485 } 486 if(ak < an) { 487 ierr = PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 488 ierr = PetscMemcpy(J_+k,aJ+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 489 k += (an-ak); 490 } 491 if(bk < bn) { 492 ierr = PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 493 ierr = PetscMemcpy(J_+k,bJ+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 494 k += (bn-bk); 495 } 496 PetscFunctionReturn(0); 497 } 498 499 500 #undef __FUNCT__ 501 #define __FUNCT__ "PetscProcessTree" 502 /*@ 503 PetscProcessTree - Prepares tree data to be displayed graphically 504 505 Not Collective 506 507 Input Parameters: 508 + n - number of values 509 . mask - indicates those entries in the tree, location 0 is always masked 510 - parentid - indicates the parent of each entry 511 512 Output Parameters: 513 + Nlevels - the number of levels 514 . Level - for each node tells its level 515 . Levelcnts - the number of nodes on each level 516 . Idbylevel - a list of ids on each of the levels, first level followed by second etc 517 - Column - for each id tells its column index 518 519 Level: intermediate 520 521 522 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 523 @*/ 524 PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column) 525 { 526 PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column; 527 PetscErrorCode ierr; 528 PetscBool done = PETSC_FALSE; 529 530 PetscFunctionBegin; 531 if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set"); 532 for (i=0; i<n; i++) { 533 if (mask[i]) continue; 534 if (parentid[i] == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent"); 535 if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked"); 536 } 537 538 539 for (i=0; i<n; i++) { 540 if (!mask[i]) nmask++; 541 } 542 543 /* determine the level in the tree of each node */ 544 ierr = PetscMalloc(n*sizeof(PetscInt),&level);CHKERRQ(ierr); 545 ierr = PetscMemzero(level,n*sizeof(PetscInt));CHKERRQ(ierr); 546 level[0] = 1; 547 while (!done) { 548 done = PETSC_TRUE; 549 for (i=0; i<n; i++) { 550 if (mask[i]) continue; 551 if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1; 552 else if (!level[i]) done = PETSC_FALSE; 553 } 554 } 555 for (i=0; i<n; i++) { 556 level[i]--; 557 nlevels = PetscMax(nlevels,level[i]); 558 } 559 560 /* count the number of nodes on each level and its max */ 561 ierr = PetscMalloc(nlevels*sizeof(PetscInt),&levelcnt);CHKERRQ(ierr); 562 ierr = PetscMemzero(levelcnt,nlevels*sizeof(PetscInt));CHKERRQ(ierr); 563 for (i=0; i<n; i++) { 564 if (mask[i]) continue; 565 levelcnt[level[i]-1]++; 566 } 567 for (i=0; i<nlevels;i++) { 568 levelmax = PetscMax(levelmax,levelcnt[i]); 569 } 570 571 /* for each level sort the ids by the parent id */ 572 ierr = PetscMalloc2(levelmax,PetscInt,&workid,levelmax,PetscInt,&workparentid);CHKERRQ(ierr); 573 ierr = PetscMalloc(nmask*sizeof(PetscInt),&idbylevel);CHKERRQ(ierr); 574 for (j=1; j<=nlevels;j++) { 575 cnt = 0; 576 for (i=0; i<n; i++) { 577 if (mask[i]) continue; 578 if (level[i] != j) continue; 579 workid[cnt] = i; 580 workparentid[cnt++] = parentid[i]; 581 } 582 /* PetscIntView(cnt,workparentid,0); 583 PetscIntView(cnt,workid,0); 584 ierr = PetscSortIntWithArray(cnt,workparentid,workid);CHKERRQ(ierr); 585 PetscIntView(cnt,workparentid,0); 586 PetscIntView(cnt,workid,0);*/ 587 ierr = PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));CHKERRQ(ierr); 588 tcnt += cnt; 589 } 590 if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes"); 591 ierr = PetscFree2(workid,workparentid);CHKERRQ(ierr); 592 593 /* for each node list its column */ 594 ierr = PetscMalloc(n*sizeof(PetscInt),&column);CHKERRQ(ierr); 595 cnt = 0; 596 for (j=0; j<nlevels; j++) { 597 for (i=0; i<levelcnt[j]; i++) { 598 column[idbylevel[cnt++]] = i; 599 } 600 } 601 602 *Nlevels = nlevels; 603 *Level = level; 604 *Levelcnt = levelcnt; 605 *Idbylevel = idbylevel; 606 *Column = column; 607 PetscFunctionReturn(0); 608 } 609