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 #undef __FUNCT__ 127 #define __FUNCT__ "PetscFindInt" 128 /*@ 129 PetscFindInt - Finds and integers in a sorted array of integers 130 131 Not Collective 132 133 Input Parameters: 134 + key - the integer to locate 135 . n - number of value in the array 136 - ii - array of integers 137 138 Output Parameter: 139 . loc - the location, or -1 if the value is not found 140 141 Level: intermediate 142 143 Concepts: sorting^ints 144 145 .seealso: PetscSortInt(), PetscSortIntWithArray(), PetscSortRemoveDupsInt() 146 @*/ 147 PetscErrorCode PetscFindInt(PetscInt key, PetscInt n, const PetscInt ii[], PetscInt *loc) 148 { 149 /* inclusive indices 150 0 <= imin when using truncate toward zero divide 151 imid = (imin+imax)/2; 152 imin unrestricted when using truncate toward minus infinity divide 153 imid = (imin+imax)>>1; or 154 imid = (int)floor((imin+imax)/2.0); */ 155 PetscInt imin = 0, imax = n-1; 156 157 PetscFunctionBegin; 158 PetscValidPointer(ii, 3); 159 PetscValidPointer(loc, 4); 160 /* continually narrow search until just one element remains */ 161 while(imin < imax) { 162 PetscInt imid = ((unsigned PetscInt) imin + (unsigned PetscInt) imax) >> 1; 163 164 /* code must guarantee the interval is reduced at each iteration 165 assert(imid < imax); */ 166 /* reduce the search */ 167 if (ii[imid] < key) { 168 imin = imid + 1; 169 } else { 170 imax = imid; 171 } 172 } 173 /* At exit of while: 174 if ii[] is empty, then imax < imin 175 otherwise imax == imin */ 176 /* deferred test for equality */ 177 if ((imax == imin) && (ii[imin] == key)) { 178 *loc = imin; 179 } else { 180 *loc = -1; 181 } 182 PetscFunctionReturn(0); 183 } 184 185 186 /* -----------------------------------------------------------------------*/ 187 #define SWAP2(a,b,c,d,t) {t=a;a=b;b=t;t=c;c=d;d=t;} 188 189 #undef __FUNCT__ 190 #define __FUNCT__ "PetscSortIntWithArray_Private" 191 /* 192 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 193 Assumes 0 origin for v, number of elements = right+1 (right is index of 194 right-most member). 195 */ 196 static PetscErrorCode PetscSortIntWithArray_Private(PetscInt *v,PetscInt *V,PetscInt right) 197 { 198 PetscErrorCode ierr; 199 PetscInt i,vl,last,tmp; 200 201 PetscFunctionBegin; 202 if (right <= 1) { 203 if (right == 1) { 204 if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp); 205 } 206 PetscFunctionReturn(0); 207 } 208 SWAP2(v[0],v[right/2],V[0],V[right/2],tmp); 209 vl = v[0]; 210 last = 0; 211 for (i=1; i<=right; i++) { 212 if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);} 213 } 214 SWAP2(v[0],v[last],V[0],V[last],tmp); 215 ierr = PetscSortIntWithArray_Private(v,V,last-1);CHKERRQ(ierr); 216 ierr = PetscSortIntWithArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 217 PetscFunctionReturn(0); 218 } 219 220 #undef __FUNCT__ 221 #define __FUNCT__ "PetscSortIntWithArray" 222 /*@ 223 PetscSortIntWithArray - Sorts an array of integers in place in increasing order; 224 changes a second array to match the sorted first array. 225 226 Not Collective 227 228 Input Parameters: 229 + n - number of values 230 . i - array of integers 231 - I - second array of integers 232 233 Level: intermediate 234 235 Concepts: sorting^ints with array 236 237 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt() 238 @*/ 239 PetscErrorCode PetscSortIntWithArray(PetscInt n,PetscInt i[],PetscInt Ii[]) 240 { 241 PetscErrorCode ierr; 242 PetscInt j,k,tmp,ik; 243 244 PetscFunctionBegin; 245 if (n<8) { 246 for (k=0; k<n; k++) { 247 ik = i[k]; 248 for (j=k+1; j<n; j++) { 249 if (ik > i[j]) { 250 SWAP2(i[k],i[j],Ii[k],Ii[j],tmp); 251 ik = i[k]; 252 } 253 } 254 } 255 } else { 256 ierr = PetscSortIntWithArray_Private(i,Ii,n-1);CHKERRQ(ierr); 257 } 258 PetscFunctionReturn(0); 259 } 260 261 262 #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;} 263 264 #undef __FUNCT__ 265 #define __FUNCT__ "PetscSortIntWithArrayPair_Private" 266 /* 267 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 268 Assumes 0 origin for v, number of elements = right+1 (right is index of 269 right-most member). 270 */ 271 static PetscErrorCode PetscSortIntWithArrayPair_Private(PetscInt *L,PetscInt *J, PetscInt *K, PetscInt right) 272 { 273 PetscErrorCode ierr; 274 PetscInt i,vl,last,tmp; 275 276 PetscFunctionBegin; 277 if (right <= 1) { 278 if (right == 1) { 279 if (L[0] > L[1]) SWAP3(L[0],L[1],J[0],J[1],K[0],K[1],tmp); 280 } 281 PetscFunctionReturn(0); 282 } 283 SWAP3(L[0],L[right/2],J[0],J[right/2],K[0],K[right/2],tmp); 284 vl = L[0]; 285 last = 0; 286 for (i=1; i<=right; i++) { 287 if (L[i] < vl) {last++; SWAP3(L[last],L[i],J[last],J[i],K[last],K[i],tmp);} 288 } 289 SWAP3(L[0],L[last],J[0],J[last],K[0],K[last],tmp); 290 ierr = PetscSortIntWithArrayPair_Private(L,J,K,last-1);CHKERRQ(ierr); 291 ierr = PetscSortIntWithArrayPair_Private(L+last+1,J+last+1,K+last+1,right-(last+1));CHKERRQ(ierr); 292 PetscFunctionReturn(0); 293 } 294 295 #undef __FUNCT__ 296 #define __FUNCT__ "PetscSortIntWithArrayPair" 297 /*@ 298 PetscSortIntWithArrayPair - Sorts an array of integers in place in increasing order; 299 changes a pair of integer arrays to match the sorted first array. 300 301 Not Collective 302 303 Input Parameters: 304 + n - number of values 305 . I - array of integers 306 . J - second array of integers (first array of the pair) 307 - K - third array of integers (second array of the pair) 308 309 Level: intermediate 310 311 Concepts: sorting^ints with array pair 312 313 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortIntWithArray() 314 @*/ 315 PetscErrorCode PetscSortIntWithArrayPair(PetscInt n,PetscInt *L,PetscInt *J, PetscInt *K) 316 { 317 PetscErrorCode ierr; 318 PetscInt j,k,tmp,ik; 319 320 PetscFunctionBegin; 321 if (n<8) { 322 for (k=0; k<n; k++) { 323 ik = L[k]; 324 for (j=k+1; j<n; j++) { 325 if (ik > L[j]) { 326 SWAP3(L[k],L[j],J[k],J[j],K[k],K[j],tmp); 327 ik = L[k]; 328 } 329 } 330 } 331 } else { 332 ierr = PetscSortIntWithArrayPair_Private(L,J,K,n-1);CHKERRQ(ierr); 333 } 334 PetscFunctionReturn(0); 335 } 336 337 #undef __FUNCT__ 338 #define __FUNCT__ "PetscSortMPIInt_Private" 339 /* 340 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 341 Assumes 0 origin for v, number of elements = right+1 (right is index of 342 right-most member). 343 */ 344 static void PetscSortMPIInt_Private(PetscMPIInt *v,PetscInt right) 345 { 346 PetscInt i,j; 347 PetscMPIInt pivot,tmp; 348 349 if (right <= 1) { 350 if (right == 1) { 351 if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 352 } 353 return; 354 } 355 i = MEDIAN(v,right); /* Choose a pivot */ 356 SWAP(v[0],v[i],tmp); /* Move it out of the way */ 357 pivot = v[0]; 358 for (i=0,j=right+1;;) { 359 while (++i < j && v[i] <= pivot); /* Scan from the left */ 360 while (v[--j] > pivot); /* Scan from the right */ 361 if (i >= j) break; 362 SWAP(v[i],v[j],tmp); 363 } 364 SWAP(v[0],v[j],tmp); /* Put pivot back in place. */ 365 PetscSortMPIInt_Private(v,j-1); 366 PetscSortMPIInt_Private(v+j+1,right-(j+1)); 367 } 368 369 #undef __FUNCT__ 370 #define __FUNCT__ "PetscSortMPIInt" 371 /*@ 372 PetscSortMPIInt - Sorts an array of MPI integers in place in increasing order. 373 374 Not Collective 375 376 Input Parameters: 377 + n - number of values 378 - i - array of integers 379 380 Level: intermediate 381 382 Concepts: sorting^ints 383 384 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 385 @*/ 386 PetscErrorCode PetscSortMPIInt(PetscInt n,PetscMPIInt i[]) 387 { 388 PetscInt j,k; 389 PetscMPIInt tmp,ik; 390 391 PetscFunctionBegin; 392 if (n<8) { 393 for (k=0; k<n; k++) { 394 ik = i[k]; 395 for (j=k+1; j<n; j++) { 396 if (ik > i[j]) { 397 SWAP(i[k],i[j],tmp); 398 ik = i[k]; 399 } 400 } 401 } 402 } else { 403 PetscSortMPIInt_Private(i,n-1); 404 } 405 PetscFunctionReturn(0); 406 } 407 408 #undef __FUNCT__ 409 #define __FUNCT__ "PetscSortRemoveDupsMPIInt" 410 /*@ 411 PetscSortRemoveDupsMPIInt - Sorts an array of MPI integers in place in increasing order removes all duplicate entries 412 413 Not Collective 414 415 Input Parameters: 416 + n - number of values 417 - ii - array of integers 418 419 Output Parameter: 420 . n - number of non-redundant values 421 422 Level: intermediate 423 424 Concepts: sorting^ints 425 426 .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt() 427 @*/ 428 PetscErrorCode PetscSortRemoveDupsMPIInt(PetscInt *n,PetscMPIInt ii[]) 429 { 430 PetscErrorCode ierr; 431 PetscInt i,s = 0,N = *n, b = 0; 432 433 PetscFunctionBegin; 434 ierr = PetscSortMPIInt(N,ii);CHKERRQ(ierr); 435 for (i=0; i<N-1; i++) { 436 if (ii[b+s+1] != ii[b]) {ii[b+1] = ii[b+s+1]; b++;} 437 else s++; 438 } 439 *n = N - s; 440 PetscFunctionReturn(0); 441 } 442 443 #undef __FUNCT__ 444 #define __FUNCT__ "PetscSortMPIIntWithArray_Private" 445 /* 446 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 447 Assumes 0 origin for v, number of elements = right+1 (right is index of 448 right-most member). 449 */ 450 static PetscErrorCode PetscSortMPIIntWithArray_Private(PetscMPIInt *v,PetscMPIInt *V,PetscMPIInt right) 451 { 452 PetscErrorCode ierr; 453 PetscMPIInt i,vl,last,tmp; 454 455 PetscFunctionBegin; 456 if (right <= 1) { 457 if (right == 1) { 458 if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp); 459 } 460 PetscFunctionReturn(0); 461 } 462 SWAP2(v[0],v[right/2],V[0],V[right/2],tmp); 463 vl = v[0]; 464 last = 0; 465 for (i=1; i<=right; i++) { 466 if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);} 467 } 468 SWAP2(v[0],v[last],V[0],V[last],tmp); 469 ierr = PetscSortMPIIntWithArray_Private(v,V,last-1);CHKERRQ(ierr); 470 ierr = PetscSortMPIIntWithArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 471 PetscFunctionReturn(0); 472 } 473 474 #undef __FUNCT__ 475 #define __FUNCT__ "PetscSortMPIIntWithArray" 476 /*@ 477 PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order; 478 changes a second array to match the sorted first array. 479 480 Not Collective 481 482 Input Parameters: 483 + n - number of values 484 . i - array of integers 485 - I - second array of integers 486 487 Level: intermediate 488 489 Concepts: sorting^ints with array 490 491 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt() 492 @*/ 493 PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt i[],PetscMPIInt Ii[]) 494 { 495 PetscErrorCode ierr; 496 PetscMPIInt j,k,tmp,ik; 497 498 PetscFunctionBegin; 499 if (n<8) { 500 for (k=0; k<n; k++) { 501 ik = i[k]; 502 for (j=k+1; j<n; j++) { 503 if (ik > i[j]) { 504 SWAP2(i[k],i[j],Ii[k],Ii[j],tmp); 505 ik = i[k]; 506 } 507 } 508 } 509 } else { 510 ierr = PetscSortMPIIntWithArray_Private(i,Ii,n-1);CHKERRQ(ierr); 511 } 512 PetscFunctionReturn(0); 513 } 514 515 /* -----------------------------------------------------------------------*/ 516 #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;} 517 518 #undef __FUNCT__ 519 #define __FUNCT__ "PetscSortIntWithScalarArray_Private" 520 /* 521 Modified from PetscSortIntWithArray_Private(). 522 */ 523 static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right) 524 { 525 PetscErrorCode ierr; 526 PetscInt i,vl,last,tmp; 527 PetscScalar stmp; 528 529 PetscFunctionBegin; 530 if (right <= 1) { 531 if (right == 1) { 532 if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp); 533 } 534 PetscFunctionReturn(0); 535 } 536 SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp); 537 vl = v[0]; 538 last = 0; 539 for (i=1; i<=right; i++) { 540 if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);} 541 } 542 SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp); 543 ierr = PetscSortIntWithScalarArray_Private(v,V,last-1);CHKERRQ(ierr); 544 ierr = PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 545 PetscFunctionReturn(0); 546 } 547 548 #undef __FUNCT__ 549 #define __FUNCT__ "PetscSortIntWithScalarArray" 550 /*@ 551 PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order; 552 changes a second SCALAR array to match the sorted first INTEGER array. 553 554 Not Collective 555 556 Input Parameters: 557 + n - number of values 558 . i - array of integers 559 - I - second array of scalars 560 561 Level: intermediate 562 563 Concepts: sorting^ints with array 564 565 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 566 @*/ 567 PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[]) 568 { 569 PetscErrorCode ierr; 570 PetscInt j,k,tmp,ik; 571 PetscScalar stmp; 572 573 PetscFunctionBegin; 574 if (n<8) { 575 for (k=0; k<n; k++) { 576 ik = i[k]; 577 for (j=k+1; j<n; j++) { 578 if (ik > i[j]) { 579 SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp); 580 ik = i[k]; 581 } 582 } 583 } 584 } else { 585 ierr = PetscSortIntWithScalarArray_Private(i,Ii,n-1);CHKERRQ(ierr); 586 } 587 PetscFunctionReturn(0); 588 } 589 590 591 592 #undef __FUNCT__ 593 #define __FUNCT__ "PetscMergeIntArrayPair" 594 /*@ 595 PetscMergeIntArrayPair - Merges two SORTED integer arrays along with an additional array of integers. 596 The additional arrays are the same length as sorted arrays and are merged 597 in the order determined by the merging of the sorted pair. 598 599 Not Collective 600 601 Input Parameters: 602 + an - number of values in the first array 603 . aI - first sorted array of integers 604 . aJ - first additional array of integers 605 . bn - number of values in the second array 606 . bI - second array of integers 607 - bJ - second additional of integers 608 609 Output Parameters: 610 + n - number of values in the merged array (== an + bn) 611 . I - merged sorted array 612 - J - merged additional array 613 614 Level: intermediate 615 616 Concepts: merging^arrays 617 618 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 619 @*/ 620 PetscErrorCode PetscMergeIntArrayPair(PetscInt an,const PetscInt *aI, const PetscInt *aJ, PetscInt bn, const PetscInt *bI, const PetscInt *bJ, PetscInt *n, PetscInt **L, PetscInt **J) 621 { 622 PetscErrorCode ierr; 623 PetscInt n_, *L_ = *L, *J_= *J, ak, bk, k; 624 625 n_ = an + bn; 626 *n = n_; 627 if (!L_) { 628 ierr = PetscMalloc(n_*sizeof(PetscInt), L); CHKERRQ(ierr); 629 L_ = *L; 630 } 631 if (!J_){ 632 ierr = PetscMalloc(n_*sizeof(PetscInt), &J_); CHKERRQ(ierr); 633 J_ = *J; 634 } 635 k = ak = bk = 0; 636 while(ak < an && bk < bn) { 637 if (aI[ak] <= bI[bk]) { 638 L_[k] = aI[ak]; 639 J_[k] = aJ[ak]; 640 ++ak; 641 ++k; 642 } 643 else { 644 L_[k] = bI[bk]; 645 J_[k] = bJ[bk]; 646 ++bk; 647 ++k; 648 } 649 } 650 if (ak < an) { 651 ierr = PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 652 ierr = PetscMemcpy(J_+k,aJ+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 653 k += (an-ak); 654 } 655 if (bk < bn) { 656 ierr = PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 657 ierr = PetscMemcpy(J_+k,bJ+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 658 k += (bn-bk); 659 } 660 PetscFunctionReturn(0); 661 } 662 663 664 #undef __FUNCT__ 665 #define __FUNCT__ "PetscProcessTree" 666 /*@ 667 PetscProcessTree - Prepares tree data to be displayed graphically 668 669 Not Collective 670 671 Input Parameters: 672 + n - number of values 673 . mask - indicates those entries in the tree, location 0 is always masked 674 - parentid - indicates the parent of each entry 675 676 Output Parameters: 677 + Nlevels - the number of levels 678 . Level - for each node tells its level 679 . Levelcnts - the number of nodes on each level 680 . Idbylevel - a list of ids on each of the levels, first level followed by second etc 681 - Column - for each id tells its column index 682 683 Level: intermediate 684 685 686 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 687 @*/ 688 PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column) 689 { 690 PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column; 691 PetscErrorCode ierr; 692 PetscBool done = PETSC_FALSE; 693 694 PetscFunctionBegin; 695 if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set"); 696 for (i=0; i<n; i++) { 697 if (mask[i]) continue; 698 if (parentid[i] == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent"); 699 if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked"); 700 } 701 702 703 for (i=0; i<n; i++) { 704 if (!mask[i]) nmask++; 705 } 706 707 /* determine the level in the tree of each node */ 708 ierr = PetscMalloc(n*sizeof(PetscInt),&level);CHKERRQ(ierr); 709 ierr = PetscMemzero(level,n*sizeof(PetscInt));CHKERRQ(ierr); 710 level[0] = 1; 711 while (!done) { 712 done = PETSC_TRUE; 713 for (i=0; i<n; i++) { 714 if (mask[i]) continue; 715 if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1; 716 else if (!level[i]) done = PETSC_FALSE; 717 } 718 } 719 for (i=0; i<n; i++) { 720 level[i]--; 721 nlevels = PetscMax(nlevels,level[i]); 722 } 723 724 /* count the number of nodes on each level and its max */ 725 ierr = PetscMalloc(nlevels*sizeof(PetscInt),&levelcnt);CHKERRQ(ierr); 726 ierr = PetscMemzero(levelcnt,nlevels*sizeof(PetscInt));CHKERRQ(ierr); 727 for (i=0; i<n; i++) { 728 if (mask[i]) continue; 729 levelcnt[level[i]-1]++; 730 } 731 for (i=0; i<nlevels;i++) { 732 levelmax = PetscMax(levelmax,levelcnt[i]); 733 } 734 735 /* for each level sort the ids by the parent id */ 736 ierr = PetscMalloc2(levelmax,PetscInt,&workid,levelmax,PetscInt,&workparentid);CHKERRQ(ierr); 737 ierr = PetscMalloc(nmask*sizeof(PetscInt),&idbylevel);CHKERRQ(ierr); 738 for (j=1; j<=nlevels;j++) { 739 cnt = 0; 740 for (i=0; i<n; i++) { 741 if (mask[i]) continue; 742 if (level[i] != j) continue; 743 workid[cnt] = i; 744 workparentid[cnt++] = parentid[i]; 745 } 746 /* PetscIntView(cnt,workparentid,0); 747 PetscIntView(cnt,workid,0); 748 ierr = PetscSortIntWithArray(cnt,workparentid,workid);CHKERRQ(ierr); 749 PetscIntView(cnt,workparentid,0); 750 PetscIntView(cnt,workid,0);*/ 751 ierr = PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));CHKERRQ(ierr); 752 tcnt += cnt; 753 } 754 if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes"); 755 ierr = PetscFree2(workid,workparentid);CHKERRQ(ierr); 756 757 /* for each node list its column */ 758 ierr = PetscMalloc(n*sizeof(PetscInt),&column);CHKERRQ(ierr); 759 cnt = 0; 760 for (j=0; j<nlevels; j++) { 761 for (i=0; i<levelcnt[j]; i++) { 762 column[idbylevel[cnt++]] = i; 763 } 764 } 765 766 *Nlevels = nlevels; 767 *Level = level; 768 *Levelcnt = levelcnt; 769 *Idbylevel = idbylevel; 770 *Column = column; 771 PetscFunctionReturn(0); 772 } 773