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