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 #undef __FUNCT__ 279 #define __FUNCT__ "PetscSortMPIInt_Private" 280 /* 281 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 282 Assumes 0 origin for v, number of elements = right+1 (right is index of 283 right-most member). 284 */ 285 static void PetscSortMPIInt_Private(PetscMPIInt *v,PetscInt right) 286 { 287 PetscInt i,j; 288 PetscMPIInt pivot,tmp; 289 290 if (right <= 1) { 291 if (right == 1) { 292 if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 293 } 294 return; 295 } 296 i = MEDIAN(v,right); /* Choose a pivot */ 297 SWAP(v[0],v[i],tmp); /* Move it out of the way */ 298 pivot = v[0]; 299 for (i=0,j=right+1;;) { 300 while (++i < j && v[i] <= pivot); /* Scan from the left */ 301 while (v[--j] > pivot); /* Scan from the right */ 302 if (i >= j) break; 303 SWAP(v[i],v[j],tmp); 304 } 305 SWAP(v[0],v[j],tmp); /* Put pivot back in place. */ 306 PetscSortMPIInt_Private(v,j-1); 307 PetscSortMPIInt_Private(v+j+1,right-(j+1)); 308 } 309 310 #undef __FUNCT__ 311 #define __FUNCT__ "PetscSortMPIInt" 312 /*@ 313 PetscSortMPIInt - Sorts an array of MPI integers in place in increasing order. 314 315 Not Collective 316 317 Input Parameters: 318 + n - number of values 319 - i - array of integers 320 321 Level: intermediate 322 323 Concepts: sorting^ints 324 325 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 326 @*/ 327 PetscErrorCode PetscSortMPIInt(PetscInt n,PetscMPIInt i[]) 328 { 329 PetscInt j,k; 330 PetscMPIInt tmp,ik; 331 332 PetscFunctionBegin; 333 if (n<8) { 334 for (k=0; k<n; k++) { 335 ik = i[k]; 336 for (j=k+1; j<n; j++) { 337 if (ik > i[j]) { 338 SWAP(i[k],i[j],tmp); 339 ik = i[k]; 340 } 341 } 342 } 343 } else { 344 PetscSortMPIInt_Private(i,n-1); 345 } 346 PetscFunctionReturn(0); 347 } 348 349 #undef __FUNCT__ 350 #define __FUNCT__ "PetscSortRemoveDupsMPIInt" 351 /*@ 352 PetscSortRemoveDupsMPIInt - Sorts an array of MPI integers in place in increasing order removes all duplicate entries 353 354 Not Collective 355 356 Input Parameters: 357 + n - number of values 358 - ii - array of integers 359 360 Output Parameter: 361 . n - number of non-redundant values 362 363 Level: intermediate 364 365 Concepts: sorting^ints 366 367 .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt() 368 @*/ 369 PetscErrorCode PetscSortRemoveDupsMPIInt(PetscInt *n,PetscMPIInt ii[]) 370 { 371 PetscErrorCode ierr; 372 PetscInt i,s = 0,N = *n, b = 0; 373 374 PetscFunctionBegin; 375 ierr = PetscSortMPIInt(N,ii);CHKERRQ(ierr); 376 for (i=0; i<N-1; i++) { 377 if (ii[b+s+1] != ii[b]) {ii[b+1] = ii[b+s+1]; b++;} 378 else s++; 379 } 380 *n = N - s; 381 PetscFunctionReturn(0); 382 } 383 384 #undef __FUNCT__ 385 #define __FUNCT__ "PetscSortMPIIntWithArray_Private" 386 /* 387 A simple version of quicksort; taken from Kernighan and Ritchie, page 87. 388 Assumes 0 origin for v, number of elements = right+1 (right is index of 389 right-most member). 390 */ 391 static PetscErrorCode PetscSortMPIIntWithArray_Private(PetscMPIInt *v,PetscMPIInt *V,PetscMPIInt right) 392 { 393 PetscErrorCode ierr; 394 PetscMPIInt i,vl,last,tmp; 395 396 PetscFunctionBegin; 397 if (right <= 1) { 398 if (right == 1) { 399 if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp); 400 } 401 PetscFunctionReturn(0); 402 } 403 SWAP2(v[0],v[right/2],V[0],V[right/2],tmp); 404 vl = v[0]; 405 last = 0; 406 for (i=1; i<=right; i++) { 407 if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);} 408 } 409 SWAP2(v[0],v[last],V[0],V[last],tmp); 410 ierr = PetscSortMPIIntWithArray_Private(v,V,last-1);CHKERRQ(ierr); 411 ierr = PetscSortMPIIntWithArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 412 PetscFunctionReturn(0); 413 } 414 415 #undef __FUNCT__ 416 #define __FUNCT__ "PetscSortMPIIntWithArray" 417 /*@ 418 PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order; 419 changes a second array to match the sorted first array. 420 421 Not Collective 422 423 Input Parameters: 424 + n - number of values 425 . i - array of integers 426 - I - second array of integers 427 428 Level: intermediate 429 430 Concepts: sorting^ints with array 431 432 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt() 433 @*/ 434 PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt i[],PetscMPIInt Ii[]) 435 { 436 PetscErrorCode ierr; 437 PetscMPIInt j,k,tmp,ik; 438 439 PetscFunctionBegin; 440 if (n<8) { 441 for (k=0; k<n; k++) { 442 ik = i[k]; 443 for (j=k+1; j<n; j++) { 444 if (ik > i[j]) { 445 SWAP2(i[k],i[j],Ii[k],Ii[j],tmp); 446 ik = i[k]; 447 } 448 } 449 } 450 } else { 451 ierr = PetscSortMPIIntWithArray_Private(i,Ii,n-1);CHKERRQ(ierr); 452 } 453 PetscFunctionReturn(0); 454 } 455 456 /* -----------------------------------------------------------------------*/ 457 #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;} 458 459 #undef __FUNCT__ 460 #define __FUNCT__ "PetscSortIntWithScalarArray_Private" 461 /* 462 Modified from PetscSortIntWithArray_Private(). 463 */ 464 static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right) 465 { 466 PetscErrorCode ierr; 467 PetscInt i,vl,last,tmp; 468 PetscScalar stmp; 469 470 PetscFunctionBegin; 471 if (right <= 1) { 472 if (right == 1) { 473 if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp); 474 } 475 PetscFunctionReturn(0); 476 } 477 SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp); 478 vl = v[0]; 479 last = 0; 480 for (i=1; i<=right; i++) { 481 if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);} 482 } 483 SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp); 484 ierr = PetscSortIntWithScalarArray_Private(v,V,last-1);CHKERRQ(ierr); 485 ierr = PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));CHKERRQ(ierr); 486 PetscFunctionReturn(0); 487 } 488 489 #undef __FUNCT__ 490 #define __FUNCT__ "PetscSortIntWithScalarArray" 491 /*@ 492 PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order; 493 changes a second SCALAR array to match the sorted first INTEGER array. 494 495 Not Collective 496 497 Input Parameters: 498 + n - number of values 499 . i - array of integers 500 - I - second array of scalars 501 502 Level: intermediate 503 504 Concepts: sorting^ints with array 505 506 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 507 @*/ 508 PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[]) 509 { 510 PetscErrorCode ierr; 511 PetscInt j,k,tmp,ik; 512 PetscScalar stmp; 513 514 PetscFunctionBegin; 515 if (n<8) { 516 for (k=0; k<n; k++) { 517 ik = i[k]; 518 for (j=k+1; j<n; j++) { 519 if (ik > i[j]) { 520 SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp); 521 ik = i[k]; 522 } 523 } 524 } 525 } else { 526 ierr = PetscSortIntWithScalarArray_Private(i,Ii,n-1);CHKERRQ(ierr); 527 } 528 PetscFunctionReturn(0); 529 } 530 531 532 533 #undef __FUNCT__ 534 #define __FUNCT__ "PetscMergeIntArrayPair" 535 /*@ 536 PetscMergeIntArrayPair - Merges two SORTED integer arrays along with an additional array of integers. 537 The additional arrays are the same length as sorted arrays and are merged 538 in the order determined by the merging of the sorted pair. 539 540 Not Collective 541 542 Input Parameters: 543 + an - number of values in the first array 544 . aI - first sorted array of integers 545 . aJ - first additional array of integers 546 . bn - number of values in the second array 547 . bI - second array of integers 548 - bJ - second additional of integers 549 550 Output Parameters: 551 + n - number of values in the merged array (== an + bn) 552 . I - merged sorted array 553 - J - merged additional array 554 555 Level: intermediate 556 557 Concepts: merging^arrays 558 559 .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray() 560 @*/ 561 PetscErrorCode PetscMergeIntArrayPair(PetscInt an,const PetscInt *aI, const PetscInt *aJ, PetscInt bn, const PetscInt *bI, const PetscInt *bJ, PetscInt *n, PetscInt **L, PetscInt **J) 562 { 563 PetscErrorCode ierr; 564 PetscInt n_, *L_ = *L, *J_= *J, ak, bk, k; 565 566 n_ = an + bn; 567 *n = n_; 568 if (!L_) { 569 ierr = PetscMalloc(n_*sizeof(PetscInt), L); CHKERRQ(ierr); 570 L_ = *L; 571 } 572 if (!J_){ 573 ierr = PetscMalloc(n_*sizeof(PetscInt), &J_); CHKERRQ(ierr); 574 J_ = *J; 575 } 576 k = ak = bk = 0; 577 while(ak < an && bk < bn) { 578 if (aI[ak] <= bI[bk]) { 579 L_[k] = aI[ak]; 580 J_[k] = aJ[ak]; 581 ++ak; 582 ++k; 583 } 584 else { 585 L_[k] = bI[bk]; 586 J_[k] = bJ[bk]; 587 ++bk; 588 ++k; 589 } 590 } 591 if (ak < an) { 592 ierr = PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 593 ierr = PetscMemcpy(J_+k,aJ+ak,(an-ak)*sizeof(PetscInt)); CHKERRQ(ierr); 594 k += (an-ak); 595 } 596 if (bk < bn) { 597 ierr = PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 598 ierr = PetscMemcpy(J_+k,bJ+bk,(bn-bk)*sizeof(PetscInt)); CHKERRQ(ierr); 599 k += (bn-bk); 600 } 601 PetscFunctionReturn(0); 602 } 603 604 605 #undef __FUNCT__ 606 #define __FUNCT__ "PetscProcessTree" 607 /*@ 608 PetscProcessTree - Prepares tree data to be displayed graphically 609 610 Not Collective 611 612 Input Parameters: 613 + n - number of values 614 . mask - indicates those entries in the tree, location 0 is always masked 615 - parentid - indicates the parent of each entry 616 617 Output Parameters: 618 + Nlevels - the number of levels 619 . Level - for each node tells its level 620 . Levelcnts - the number of nodes on each level 621 . Idbylevel - a list of ids on each of the levels, first level followed by second etc 622 - Column - for each id tells its column index 623 624 Level: intermediate 625 626 627 .seealso: PetscSortReal(), PetscSortIntWithPermutation() 628 @*/ 629 PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column) 630 { 631 PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column; 632 PetscErrorCode ierr; 633 PetscBool done = PETSC_FALSE; 634 635 PetscFunctionBegin; 636 if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set"); 637 for (i=0; i<n; i++) { 638 if (mask[i]) continue; 639 if (parentid[i] == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent"); 640 if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked"); 641 } 642 643 644 for (i=0; i<n; i++) { 645 if (!mask[i]) nmask++; 646 } 647 648 /* determine the level in the tree of each node */ 649 ierr = PetscMalloc(n*sizeof(PetscInt),&level);CHKERRQ(ierr); 650 ierr = PetscMemzero(level,n*sizeof(PetscInt));CHKERRQ(ierr); 651 level[0] = 1; 652 while (!done) { 653 done = PETSC_TRUE; 654 for (i=0; i<n; i++) { 655 if (mask[i]) continue; 656 if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1; 657 else if (!level[i]) done = PETSC_FALSE; 658 } 659 } 660 for (i=0; i<n; i++) { 661 level[i]--; 662 nlevels = PetscMax(nlevels,level[i]); 663 } 664 665 /* count the number of nodes on each level and its max */ 666 ierr = PetscMalloc(nlevels*sizeof(PetscInt),&levelcnt);CHKERRQ(ierr); 667 ierr = PetscMemzero(levelcnt,nlevels*sizeof(PetscInt));CHKERRQ(ierr); 668 for (i=0; i<n; i++) { 669 if (mask[i]) continue; 670 levelcnt[level[i]-1]++; 671 } 672 for (i=0; i<nlevels;i++) { 673 levelmax = PetscMax(levelmax,levelcnt[i]); 674 } 675 676 /* for each level sort the ids by the parent id */ 677 ierr = PetscMalloc2(levelmax,PetscInt,&workid,levelmax,PetscInt,&workparentid);CHKERRQ(ierr); 678 ierr = PetscMalloc(nmask*sizeof(PetscInt),&idbylevel);CHKERRQ(ierr); 679 for (j=1; j<=nlevels;j++) { 680 cnt = 0; 681 for (i=0; i<n; i++) { 682 if (mask[i]) continue; 683 if (level[i] != j) continue; 684 workid[cnt] = i; 685 workparentid[cnt++] = parentid[i]; 686 } 687 /* PetscIntView(cnt,workparentid,0); 688 PetscIntView(cnt,workid,0); 689 ierr = PetscSortIntWithArray(cnt,workparentid,workid);CHKERRQ(ierr); 690 PetscIntView(cnt,workparentid,0); 691 PetscIntView(cnt,workid,0);*/ 692 ierr = PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));CHKERRQ(ierr); 693 tcnt += cnt; 694 } 695 if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes"); 696 ierr = PetscFree2(workid,workparentid);CHKERRQ(ierr); 697 698 /* for each node list its column */ 699 ierr = PetscMalloc(n*sizeof(PetscInt),&column);CHKERRQ(ierr); 700 cnt = 0; 701 for (j=0; j<nlevels; j++) { 702 for (i=0; i<levelcnt[j]; i++) { 703 column[idbylevel[cnt++]] = i; 704 } 705 } 706 707 *Nlevels = nlevels; 708 *Level = level; 709 *Levelcnt = levelcnt; 710 *Idbylevel = idbylevel; 711 *Column = column; 712 PetscFunctionReturn(0); 713 } 714