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