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