1 /**********************************ivec.c************************************** 2 3 Author: Henry M. Tufo III 4 5 e-mail: hmt@cs.brown.edu 6 7 snail-mail: 8 Division of Applied Mathematics 9 Brown University 10 Providence, RI 02912 11 12 Last Modification: 13 6.21.97 14 ***********************************ivec.c*************************************/ 15 16 #include <../src/ksp/pc/impls/tfs/tfs.h> 17 18 /* sorting args ivec.c ivec.c ... */ 19 #define SORT_OPT 6 20 #define SORT_STACK 50000 21 22 /* allocate an address and size stack for sorter(s) */ 23 static void *offset_stack[2 * SORT_STACK]; 24 static PetscInt size_stack[SORT_STACK]; 25 26 /***********************************ivec.c*************************************/ 27 PetscInt *PCTFS_ivec_copy(PetscInt *arg1, PetscInt *arg2, PetscInt n) 28 { 29 while (n--) *arg1++ = *arg2++; 30 return arg1; 31 } 32 33 /***********************************ivec.c*************************************/ 34 PetscErrorCode PCTFS_ivec_zero(PetscInt *arg1, PetscInt n) 35 { 36 PetscFunctionBegin; 37 while (n--) *arg1++ = 0; 38 PetscFunctionReturn(PETSC_SUCCESS); 39 } 40 41 /***********************************ivec.c*************************************/ 42 PetscErrorCode PCTFS_ivec_set(PetscInt *arg1, PetscInt arg2, PetscInt n) 43 { 44 PetscFunctionBegin; 45 while (n--) *arg1++ = arg2; 46 PetscFunctionReturn(PETSC_SUCCESS); 47 } 48 49 /***********************************ivec.c*************************************/ 50 PetscErrorCode PCTFS_ivec_max(PetscInt *arg1, PetscInt *arg2, PetscInt n) 51 { 52 PetscFunctionBegin; 53 while (n--) { 54 *arg1 = PetscMax(*arg1, *arg2); 55 arg1++; 56 arg2++; 57 } 58 PetscFunctionReturn(PETSC_SUCCESS); 59 } 60 61 /***********************************ivec.c*************************************/ 62 PetscErrorCode PCTFS_ivec_min(PetscInt *arg1, PetscInt *arg2, PetscInt n) 63 { 64 PetscFunctionBegin; 65 while (n--) { 66 *(arg1) = PetscMin(*arg1, *arg2); 67 arg1++; 68 arg2++; 69 } 70 PetscFunctionReturn(PETSC_SUCCESS); 71 } 72 73 /***********************************ivec.c*************************************/ 74 PetscErrorCode PCTFS_ivec_mult(PetscInt *arg1, PetscInt *arg2, PetscInt n) 75 { 76 PetscFunctionBegin; 77 while (n--) *arg1++ *= *arg2++; 78 PetscFunctionReturn(PETSC_SUCCESS); 79 } 80 81 /***********************************ivec.c*************************************/ 82 PetscErrorCode PCTFS_ivec_add(PetscInt *arg1, PetscInt *arg2, PetscInt n) 83 { 84 PetscFunctionBegin; 85 while (n--) *arg1++ += *arg2++; 86 PetscFunctionReturn(PETSC_SUCCESS); 87 } 88 89 /***********************************ivec.c*************************************/ 90 PetscErrorCode PCTFS_ivec_lxor(PetscInt *arg1, PetscInt *arg2, PetscInt n) 91 { 92 PetscFunctionBegin; 93 while (n--) { 94 *arg1 = (*arg1 || *arg2) && !(*arg1 && *arg2); 95 arg1++; 96 arg2++; 97 } 98 PetscFunctionReturn(PETSC_SUCCESS); 99 } 100 101 /***********************************ivec.c*************************************/ 102 PetscErrorCode PCTFS_ivec_xor(PetscInt *arg1, PetscInt *arg2, PetscInt n) 103 { 104 PetscFunctionBegin; 105 while (n--) *arg1++ ^= *arg2++; 106 PetscFunctionReturn(PETSC_SUCCESS); 107 } 108 109 /***********************************ivec.c*************************************/ 110 PetscErrorCode PCTFS_ivec_or(PetscInt *arg1, PetscInt *arg2, PetscInt n) 111 { 112 PetscFunctionBegin; 113 while (n--) *arg1++ |= *arg2++; 114 PetscFunctionReturn(PETSC_SUCCESS); 115 } 116 117 /***********************************ivec.c*************************************/ 118 PetscErrorCode PCTFS_ivec_lor(PetscInt *arg1, PetscInt *arg2, PetscInt n) 119 { 120 PetscFunctionBegin; 121 while (n--) { 122 *arg1 = (*arg1 || *arg2); 123 arg1++; 124 arg2++; 125 } 126 PetscFunctionReturn(PETSC_SUCCESS); 127 } 128 129 /***********************************ivec.c*************************************/ 130 PetscErrorCode PCTFS_ivec_and(PetscInt *arg1, PetscInt *arg2, PetscInt n) 131 { 132 PetscFunctionBegin; 133 while (n--) *arg1++ &= *arg2++; 134 PetscFunctionReturn(PETSC_SUCCESS); 135 } 136 137 /***********************************ivec.c*************************************/ 138 PetscErrorCode PCTFS_ivec_land(PetscInt *arg1, PetscInt *arg2, PetscInt n) 139 { 140 PetscFunctionBegin; 141 while (n--) { 142 *arg1 = (*arg1 && *arg2); 143 arg1++; 144 arg2++; 145 } 146 PetscFunctionReturn(PETSC_SUCCESS); 147 } 148 149 /***********************************ivec.c*************************************/ 150 PetscErrorCode PCTFS_ivec_and3(PetscInt *arg1, PetscInt *arg2, PetscInt *arg3, PetscInt n) 151 { 152 PetscFunctionBegin; 153 while (n--) *arg1++ = (*arg2++ & *arg3++); 154 PetscFunctionReturn(PETSC_SUCCESS); 155 } 156 157 /***********************************ivec.c*************************************/ 158 PetscInt PCTFS_ivec_sum(PetscInt *arg1, PetscInt n) 159 { 160 PetscInt tmp = 0; 161 while (n--) tmp += *arg1++; 162 return tmp; 163 } 164 165 /***********************************ivec.c*************************************/ 166 PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *arg1, PetscInt *arg2, PetscInt n, ...) 167 { 168 PetscInt i, j, type; 169 PetscInt *arg3; 170 va_list ap; 171 172 PetscFunctionBegin; 173 va_start(ap, n); 174 arg3 = va_arg(ap, PetscInt *); 175 va_end(ap); 176 177 /* LATER: if we're really motivated we can sort and then unsort */ 178 for (i = 0; i < n;) { 179 /* clump 'em for now */ 180 j = i + 1; 181 type = arg3[i]; 182 while ((j < n) && (arg3[j] == type)) j++; 183 184 /* how many together */ 185 j -= i; 186 187 /* call appropriate ivec function */ 188 if (type == GL_MAX) PetscCall(PCTFS_ivec_max(arg1, arg2, j)); 189 else if (type == GL_MIN) PetscCall(PCTFS_ivec_min(arg1, arg2, j)); 190 else if (type == GL_MULT) PetscCall(PCTFS_ivec_mult(arg1, arg2, j)); 191 else if (type == GL_ADD) PetscCall(PCTFS_ivec_add(arg1, arg2, j)); 192 else if (type == GL_B_XOR) PetscCall(PCTFS_ivec_xor(arg1, arg2, j)); 193 else if (type == GL_B_OR) PetscCall(PCTFS_ivec_or(arg1, arg2, j)); 194 else if (type == GL_B_AND) PetscCall(PCTFS_ivec_and(arg1, arg2, j)); 195 else if (type == GL_L_XOR) PetscCall(PCTFS_ivec_lxor(arg1, arg2, j)); 196 else if (type == GL_L_OR) PetscCall(PCTFS_ivec_lor(arg1, arg2, j)); 197 else if (type == GL_L_AND) PetscCall(PCTFS_ivec_land(arg1, arg2, j)); 198 else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "unrecognized type passed to PCTFS_ivec_non_uniform()!"); 199 200 arg1 += j; 201 arg2 += j; 202 i += j; 203 } 204 PetscFunctionReturn(PETSC_SUCCESS); 205 } 206 207 /***********************************ivec.c*************************************/ 208 vfp PCTFS_ivec_fct_addr(PetscInt type) 209 { 210 if (type == NON_UNIFORM) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_non_uniform; 211 else if (type == GL_MAX) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_max; 212 else if (type == GL_MIN) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_min; 213 else if (type == GL_MULT) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_mult; 214 else if (type == GL_ADD) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_add; 215 else if (type == GL_B_XOR) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_xor; 216 else if (type == GL_B_OR) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_or; 217 else if (type == GL_B_AND) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_and; 218 else if (type == GL_L_XOR) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_lxor; 219 else if (type == GL_L_OR) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_lor; 220 else if (type == GL_L_AND) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_ivec_land; 221 222 /* catch all ... not good if we get here */ 223 return NULL; 224 } 225 226 /******************************************************************************/ 227 PetscErrorCode PCTFS_ivec_sort(PetscInt *ar, PetscInt size) 228 { 229 PetscInt *pi, *pj, temp; 230 PetscInt **top_a = (PetscInt **)offset_stack; 231 PetscInt *top_s = size_stack, *bottom_s = size_stack; 232 233 PetscFunctionBegin; 234 /* we're really interested in the offset of the last element */ 235 /* ==> length of the list is now size + 1 */ 236 size--; 237 238 /* do until we're done ... return when stack is exhausted */ 239 for (;;) { 240 /* if list is large enough use quicksort partition exchange code */ 241 if (size > SORT_OPT) { 242 /* start up pointer at element 1 and down at size */ 243 pi = ar + 1; 244 pj = ar + size; 245 246 /* find middle element in list and swap w/ element 1 */ 247 SWAP(*(ar + (size >> 1)), *pi); 248 249 /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */ 250 /* note ==> pivot_value in index 0 */ 251 if (*pi > *pj) SWAP(*pi, *pj); 252 if (*ar > *pj) SWAP(*ar, *pj); 253 else if (*pi > *ar) SWAP(*(ar), *(ar + 1)); 254 255 /* partition about pivot_value ... */ 256 /* note lists of length 2 are not guaranteed to be sorted */ 257 for (;;) { 258 /* walk up ... and down ... swap if equal to pivot! */ 259 do pi++; 260 while (*pi < *ar); 261 do pj--; 262 while (*pj > *ar); 263 264 /* if we've crossed we're done */ 265 if (pj < pi) break; 266 267 /* else swap */ 268 SWAP(*pi, *pj); 269 } 270 271 /* place pivot_value in it's correct location */ 272 SWAP(*ar, *pj); 273 274 /* test stack_size to see if we've exhausted our stack */ 275 PetscCheck(top_s - bottom_s < SORT_STACK, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_ivec_sort() :: STACK EXHAUSTED!!!"); 276 277 /* push right hand child iff length > 1 */ 278 if ((*top_s = size - ((PetscInt)(pi - ar)))) { 279 *(top_a++) = pi; 280 size -= *top_s + 2; 281 top_s++; 282 } else if (size -= *top_s + 2) 283 ; /* set up for next loop iff there is something to do */ 284 else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */ ar = *(--top_a); 285 size = *(--top_s); 286 } 287 } else { /* else sort small list directly then pop another off stack */ 288 289 /* insertion sort for bottom */ 290 for (pj = ar + 1; pj <= ar + size; pj++) { 291 temp = *pj; 292 for (pi = pj - 1; pi >= ar; pi--) { 293 if (*pi <= temp) break; 294 *(pi + 1) = *pi; 295 } 296 *(pi + 1) = temp; 297 } 298 299 /* check to see if stack is exhausted ==> DONE */ 300 if (top_s == bottom_s) PetscFunctionReturn(PETSC_SUCCESS); 301 302 /* else pop another list from the stack */ 303 ar = *(--top_a); 304 size = *(--top_s); 305 } 306 } 307 } 308 309 /******************************************************************************/ 310 PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *ar, PetscInt *ar2, PetscInt size) 311 { 312 PetscInt *pi, *pj, temp, temp2; 313 PetscInt **top_a = (PetscInt **)offset_stack; 314 PetscInt *top_s = size_stack, *bottom_s = size_stack; 315 PetscInt *pi2, *pj2; 316 PetscInt mid; 317 318 PetscFunctionBegin; 319 /* we're really interested in the offset of the last element */ 320 /* ==> length of the list is now size + 1 */ 321 size--; 322 323 /* do until we're done ... return when stack is exhausted */ 324 for (;;) { 325 /* if list is large enough use quicksort partition exchange code */ 326 if (size > SORT_OPT) { 327 /* start up pointer at element 1 and down at size */ 328 mid = size >> 1; 329 pi = ar + 1; 330 pj = ar + mid; 331 pi2 = ar2 + 1; 332 pj2 = ar2 + mid; 333 334 /* find middle element in list and swap w/ element 1 */ 335 SWAP(*pi, *pj); 336 SWAP(*pi2, *pj2); 337 338 /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */ 339 /* note ==> pivot_value in index 0 */ 340 pj = ar + size; 341 pj2 = ar2 + size; 342 if (*pi > *pj) { 343 SWAP(*pi, *pj); 344 SWAP(*pi2, *pj2); 345 } 346 if (*ar > *pj) { 347 SWAP(*ar, *pj); 348 SWAP(*ar2, *pj2); 349 } else if (*pi > *ar) { 350 SWAP(*(ar), *(ar + 1)); 351 SWAP(*(ar2), *(ar2 + 1)); 352 } 353 354 /* partition about pivot_value ... */ 355 /* note lists of length 2 are not guaranteed to be sorted */ 356 for (;;) { 357 /* walk up ... and down ... swap if equal to pivot! */ 358 do { 359 pi++; 360 pi2++; 361 } while (*pi < *ar); 362 do { 363 pj--; 364 pj2--; 365 } while (*pj > *ar); 366 367 /* if we've crossed we're done */ 368 if (pj < pi) break; 369 370 /* else swap */ 371 SWAP(*pi, *pj); 372 SWAP(*pi2, *pj2); 373 } 374 375 /* place pivot_value in it's correct location */ 376 SWAP(*ar, *pj); 377 SWAP(*ar2, *pj2); 378 379 /* test stack_size to see if we've exhausted our stack */ 380 PetscCheck(top_s - bottom_s < SORT_STACK, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_ivec_sort_companion() :: STACK EXHAUSTED!!!"); 381 382 /* push right hand child iff length > 1 */ 383 if ((*top_s = size - ((PetscInt)(pi - ar)))) { 384 *(top_a++) = pi; 385 *(top_a++) = pi2; 386 size -= *top_s + 2; 387 top_s++; 388 } else if (size -= *top_s + 2) 389 ; /* set up for next loop iff there is something to do */ 390 else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */ ar2 = *(--top_a); 391 ar = *(--top_a); 392 size = *(--top_s); 393 } 394 } else { /* else sort small list directly then pop another off stack */ 395 396 /* insertion sort for bottom */ 397 for (pj = ar + 1, pj2 = ar2 + 1; pj <= ar + size; pj++, pj2++) { 398 temp = *pj; 399 temp2 = *pj2; 400 for (pi = pj - 1, pi2 = pj2 - 1; pi >= ar; pi--, pi2--) { 401 if (*pi <= temp) break; 402 *(pi + 1) = *pi; 403 *(pi2 + 1) = *pi2; 404 } 405 *(pi + 1) = temp; 406 *(pi2 + 1) = temp2; 407 } 408 409 /* check to see if stack is exhausted ==> DONE */ 410 if (top_s == bottom_s) PetscFunctionReturn(PETSC_SUCCESS); 411 412 /* else pop another list from the stack */ 413 ar2 = *(--top_a); 414 ar = *(--top_a); 415 size = *(--top_s); 416 } 417 } 418 } 419 420 /******************************************************************************/ 421 PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *ar, PetscInt **ar2, PetscInt size) 422 { 423 PetscInt *pi, *pj, temp, *ptr; 424 PetscInt **top_a = (PetscInt **)offset_stack; 425 PetscInt *top_s = size_stack, *bottom_s = size_stack; 426 PetscInt **pi2, **pj2; 427 PetscInt mid; 428 429 PetscFunctionBegin; 430 /* we're really interested in the offset of the last element */ 431 /* ==> length of the list is now size + 1 */ 432 size--; 433 434 /* do until we're done ... return when stack is exhausted */ 435 for (;;) { 436 /* if list is large enough use quicksort partition exchange code */ 437 if (size > SORT_OPT) { 438 /* start up pointer at element 1 and down at size */ 439 mid = size >> 1; 440 pi = ar + 1; 441 pj = ar + mid; 442 pi2 = ar2 + 1; 443 pj2 = ar2 + mid; 444 445 /* find middle element in list and swap w/ element 1 */ 446 SWAP(*pi, *pj); 447 P_SWAP(*pi2, *pj2); 448 449 /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */ 450 /* note ==> pivot_value in index 0 */ 451 pj = ar + size; 452 pj2 = ar2 + size; 453 if (*pi > *pj) { 454 SWAP(*pi, *pj); 455 P_SWAP(*pi2, *pj2); 456 } 457 if (*ar > *pj) { 458 SWAP(*ar, *pj); 459 P_SWAP(*ar2, *pj2); 460 } else if (*pi > *ar) { 461 SWAP(*(ar), *(ar + 1)); 462 P_SWAP(*(ar2), *(ar2 + 1)); 463 } 464 465 /* partition about pivot_value ... */ 466 /* note lists of length 2 are not guaranteed to be sorted */ 467 for (;;) { 468 /* walk up ... and down ... swap if equal to pivot! */ 469 do { 470 pi++; 471 pi2++; 472 } while (*pi < *ar); 473 do { 474 pj--; 475 pj2--; 476 } while (*pj > *ar); 477 478 /* if we've crossed we're done */ 479 if (pj < pi) break; 480 481 /* else swap */ 482 SWAP(*pi, *pj); 483 P_SWAP(*pi2, *pj2); 484 } 485 486 /* place pivot_value in it's correct location */ 487 SWAP(*ar, *pj); 488 P_SWAP(*ar2, *pj2); 489 490 /* test stack_size to see if we've exhausted our stack */ 491 PetscCheck(top_s - bottom_s < SORT_STACK, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_ivec_sort_companion_hack() :: STACK EXHAUSTED!!!"); 492 493 /* push right hand child iff length > 1 */ 494 if ((*top_s = size - ((PetscInt)(pi - ar)))) { 495 *(top_a++) = pi; 496 *(top_a++) = (PetscInt *)pi2; 497 size -= *top_s + 2; 498 top_s++; 499 } else if (size -= *top_s + 2) 500 ; /* set up for next loop iff there is something to do */ 501 else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */ ar2 = (PetscInt **)*(--top_a); 502 ar = *(--top_a); 503 size = *(--top_s); 504 } 505 } else { /* else sort small list directly then pop another off stack */ 506 /* insertion sort for bottom */ 507 for (pj = ar + 1, pj2 = ar2 + 1; pj <= ar + size; pj++, pj2++) { 508 temp = *pj; 509 ptr = *pj2; 510 for (pi = pj - 1, pi2 = pj2 - 1; pi >= ar; pi--, pi2--) { 511 if (*pi <= temp) break; 512 *(pi + 1) = *pi; 513 *(pi2 + 1) = *pi2; 514 } 515 *(pi + 1) = temp; 516 *(pi2 + 1) = ptr; 517 } 518 519 /* check to see if stack is exhausted ==> DONE */ 520 if (top_s == bottom_s) PetscFunctionReturn(PETSC_SUCCESS); 521 522 /* else pop another list from the stack */ 523 ar2 = (PetscInt **)*(--top_a); 524 ar = *(--top_a); 525 size = *(--top_s); 526 } 527 } 528 } 529 530 /******************************************************************************/ 531 PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type) 532 { 533 PetscFunctionBegin; 534 if (type == SORT_INTEGER) { 535 if (ar2) PetscCall(PCTFS_ivec_sort_companion((PetscInt *)ar1, (PetscInt *)ar2, size)); 536 else PetscCall(PCTFS_ivec_sort((PetscInt *)ar1, size)); 537 } else if (type == SORT_INT_PTR) { 538 if (ar2) PetscCall(PCTFS_ivec_sort_companion_hack((PetscInt *)ar1, (PetscInt **)ar2, size)); 539 else PetscCall(PCTFS_ivec_sort((PetscInt *)ar1, size)); 540 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_SMI_sort only does SORT_INTEGER!"); 541 PetscFunctionReturn(PETSC_SUCCESS); 542 } 543 544 /***********************************ivec.c*************************************/ 545 PetscInt PCTFS_ivec_linear_search(PetscInt item, PetscInt *list, PetscInt n) 546 { 547 PetscInt tmp = n - 1; 548 549 while (n--) { 550 if (*list++ == item) return tmp - n; 551 } 552 return -1; 553 } 554 555 /***********************************ivec.c*************************************/ 556 PetscInt PCTFS_ivec_binary_search(PetscInt item, PetscInt *list, PetscInt rh) 557 { 558 PetscInt mid, lh = 0; 559 560 rh--; 561 while (lh <= rh) { 562 mid = (lh + rh) >> 1; 563 if (*(list + mid) == item) return mid; 564 if (*(list + mid) > item) rh = mid - 1; 565 else lh = mid + 1; 566 } 567 return -1; 568 } 569 570 /*********************************ivec.c*************************************/ 571 PetscErrorCode PCTFS_rvec_copy(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 572 { 573 PetscFunctionBegin; 574 while (n--) *arg1++ = *arg2++; 575 PetscFunctionReturn(PETSC_SUCCESS); 576 } 577 578 /*********************************ivec.c*************************************/ 579 PetscErrorCode PCTFS_rvec_zero(PetscScalar *arg1, PetscInt n) 580 { 581 PetscFunctionBegin; 582 while (n--) *arg1++ = 0.0; 583 PetscFunctionReturn(PETSC_SUCCESS); 584 } 585 586 /***********************************ivec.c*************************************/ 587 PetscErrorCode PCTFS_rvec_one(PetscScalar *arg1, PetscInt n) 588 { 589 PetscFunctionBegin; 590 while (n--) *arg1++ = 1.0; 591 PetscFunctionReturn(PETSC_SUCCESS); 592 } 593 594 /***********************************ivec.c*************************************/ 595 PetscErrorCode PCTFS_rvec_set(PetscScalar *arg1, PetscScalar arg2, PetscInt n) 596 { 597 PetscFunctionBegin; 598 while (n--) *arg1++ = arg2; 599 PetscFunctionReturn(PETSC_SUCCESS); 600 } 601 602 /***********************************ivec.c*************************************/ 603 PetscErrorCode PCTFS_rvec_scale(PetscScalar *arg1, PetscScalar arg2, PetscInt n) 604 { 605 PetscFunctionBegin; 606 while (n--) *arg1++ *= arg2; 607 PetscFunctionReturn(PETSC_SUCCESS); 608 } 609 610 /*********************************ivec.c*************************************/ 611 PetscErrorCode PCTFS_rvec_add(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 612 { 613 PetscFunctionBegin; 614 while (n--) *arg1++ += *arg2++; 615 PetscFunctionReturn(PETSC_SUCCESS); 616 } 617 618 /*********************************ivec.c*************************************/ 619 PetscErrorCode PCTFS_rvec_mult(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 620 { 621 PetscFunctionBegin; 622 while (n--) *arg1++ *= *arg2++; 623 PetscFunctionReturn(PETSC_SUCCESS); 624 } 625 626 /*********************************ivec.c*************************************/ 627 PetscErrorCode PCTFS_rvec_max(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 628 { 629 PetscFunctionBegin; 630 while (n--) { 631 *arg1 = PetscMax(*arg1, *arg2); 632 arg1++; 633 arg2++; 634 } 635 PetscFunctionReturn(PETSC_SUCCESS); 636 } 637 638 /*********************************ivec.c*************************************/ 639 PetscErrorCode PCTFS_rvec_max_abs(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 640 { 641 PetscFunctionBegin; 642 while (n--) { 643 *arg1 = MAX_FABS(*arg1, *arg2); 644 arg1++; 645 arg2++; 646 } 647 PetscFunctionReturn(PETSC_SUCCESS); 648 } 649 650 /*********************************ivec.c*************************************/ 651 PetscErrorCode PCTFS_rvec_min(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 652 { 653 PetscFunctionBegin; 654 while (n--) { 655 *arg1 = PetscMin(*arg1, *arg2); 656 arg1++; 657 arg2++; 658 } 659 PetscFunctionReturn(PETSC_SUCCESS); 660 } 661 662 /*********************************ivec.c*************************************/ 663 PetscErrorCode PCTFS_rvec_min_abs(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 664 { 665 PetscFunctionBegin; 666 while (n--) { 667 *arg1 = MIN_FABS(*arg1, *arg2); 668 arg1++; 669 arg2++; 670 } 671 PetscFunctionReturn(PETSC_SUCCESS); 672 } 673 674 /*********************************ivec.c*************************************/ 675 static PetscErrorCode PCTFS_rvec_exists(PetscScalar *arg1, PetscScalar *arg2, PetscInt n) 676 { 677 PetscFunctionBegin; 678 while (n--) { 679 *arg1 = EXISTS(*arg1, *arg2); 680 arg1++; 681 arg2++; 682 } 683 PetscFunctionReturn(PETSC_SUCCESS); 684 } 685 686 /***********************************ivec.c*************************************/ 687 static PetscErrorCode PCTFS_rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2, PetscInt n, PetscInt *arg3) 688 { 689 PetscInt i, j, type; 690 691 PetscFunctionBegin; 692 /* LATER: if we're really motivated we can sort and then unsort */ 693 for (i = 0; i < n;) { 694 /* clump 'em for now */ 695 j = i + 1; 696 type = arg3[i]; 697 while ((j < n) && (arg3[j] == type)) j++; 698 699 /* how many together */ 700 j -= i; 701 702 /* call appropriate ivec function */ 703 if (type == GL_MAX) PetscCall(PCTFS_rvec_max(arg1, arg2, j)); 704 else if (type == GL_MIN) PetscCall(PCTFS_rvec_min(arg1, arg2, j)); 705 else if (type == GL_MULT) PetscCall(PCTFS_rvec_mult(arg1, arg2, j)); 706 else if (type == GL_ADD) PetscCall(PCTFS_rvec_add(arg1, arg2, j)); 707 else if (type == GL_MAX_ABS) PetscCall(PCTFS_rvec_max_abs(arg1, arg2, j)); 708 else if (type == GL_MIN_ABS) PetscCall(PCTFS_rvec_min_abs(arg1, arg2, j)); 709 else if (type == GL_EXISTS) PetscCall(PCTFS_rvec_exists(arg1, arg2, j)); 710 else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "unrecognized type passed to PCTFS_rvec_non_uniform()!"); 711 712 arg1 += j; 713 arg2 += j; 714 i += j; 715 } 716 PetscFunctionReturn(PETSC_SUCCESS); 717 } 718 719 /***********************************ivec.c*************************************/ 720 vfp PCTFS_rvec_fct_addr(PetscInt type) 721 { 722 if (type == NON_UNIFORM) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_non_uniform; 723 else if (type == GL_MAX) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_max; 724 else if (type == GL_MIN) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_min; 725 else if (type == GL_MULT) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_mult; 726 else if (type == GL_ADD) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_add; 727 else if (type == GL_MAX_ABS) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_max_abs; 728 else if (type == GL_MIN_ABS) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_min_abs; 729 else if (type == GL_EXISTS) return (PetscErrorCode (*)(void *, void *, PetscInt, ...))&PCTFS_rvec_exists; 730 731 /* catch all ... not good if we get here */ 732 return NULL; 733 } 734