Actual source code: ivec.c


  2: /**********************************ivec.c**************************************

  4: Author: Henry M. Tufo III

  6: e-mail: hmt@cs.brown.edu

  8: snail-mail:
  9: Division of Applied Mathematics
 10: Brown University
 11: Providence, RI 02912

 13: Last Modification:
 14: 6.21.97
 15: ***********************************ivec.c*************************************/

 17: #include <../src/ksp/pc/impls/tfs/tfs.h>

 19: /* sorting args ivec.c ivec.c ... */
 20: #define   SORT_OPT     6
 21: #define   SORT_STACK   50000

 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];

 27: /***********************************ivec.c*************************************/
 28: PetscInt *PCTFS_ivec_copy(PetscInt *arg1, PetscInt *arg2, PetscInt n)
 29: {
 30:   while (n--) *arg1++ = *arg2++;
 31:   return(arg1);
 32: }

 34: /***********************************ivec.c*************************************/
 35: PetscErrorCode PCTFS_ivec_zero(PetscInt *arg1, PetscInt n)
 36: {
 37:   while (n--) *arg1++ = 0;
 38:   return 0;
 39: }

 41: /***********************************ivec.c*************************************/
 42: PetscErrorCode PCTFS_ivec_set(PetscInt *arg1, PetscInt arg2, PetscInt n)
 43: {
 44:   while (n--) *arg1++ = arg2;
 45:   return 0;
 46: }

 48: /***********************************ivec.c*************************************/
 49: PetscErrorCode PCTFS_ivec_max(PetscInt *arg1, PetscInt *arg2, PetscInt n)
 50: {
 51:   while (n--) { *arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++; }
 52:   return 0;
 53: }

 55: /***********************************ivec.c*************************************/
 56: PetscErrorCode PCTFS_ivec_min(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 57: {
 58:   while (n--) {
 59:     *(arg1) = PetscMin(*arg1,*arg2);
 60:     arg1++;
 61:     arg2++;
 62:   }
 63:   return 0;
 64: }

 66: /***********************************ivec.c*************************************/
 67: PetscErrorCode PCTFS_ivec_mult(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 68: {
 69:   while (n--) *arg1++ *= *arg2++;
 70:   return 0;
 71: }

 73: /***********************************ivec.c*************************************/
 74: PetscErrorCode PCTFS_ivec_add(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 75: {
 76:   while (n--) *arg1++ += *arg2++;
 77:   return 0;
 78: }

 80: /***********************************ivec.c*************************************/
 81: PetscErrorCode PCTFS_ivec_lxor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 82: {
 83:   while (n--) {
 84:     *arg1=((*arg1 || *arg2) && !(*arg1 && *arg2));
 85:     arg1++;
 86:     arg2++;
 87:   }
 88:   return 0;
 89: }

 91: /***********************************ivec.c*************************************/
 92: PetscErrorCode PCTFS_ivec_xor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 93: {
 94:   while (n--) *arg1++ ^= *arg2++;
 95:   return 0;
 96: }

 98: /***********************************ivec.c*************************************/
 99: PetscErrorCode PCTFS_ivec_or(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
100: {
101:   while (n--) *arg1++ |= *arg2++;
102:   return 0;
103: }

105: /***********************************ivec.c*************************************/
106: PetscErrorCode PCTFS_ivec_lor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
107: {
108:   while (n--) {
109:     *arg1 = (*arg1 || *arg2);
110:     arg1++;
111:     arg2++;
112:   }
113:   return 0;
114: }

116: /***********************************ivec.c*************************************/
117: PetscErrorCode PCTFS_ivec_and(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
118: {
119:   while (n--) *arg1++ &= *arg2++;
120:   return 0;
121: }

123: /***********************************ivec.c*************************************/
124: PetscErrorCode PCTFS_ivec_land(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
125: {
126:   while (n--) {
127:     *arg1 = (*arg1 && *arg2);
128:     arg1++;
129:     arg2++;
130:   }
131:   return 0;
132: }

134: /***********************************ivec.c*************************************/
135: PetscErrorCode PCTFS_ivec_and3(PetscInt *arg1,  PetscInt *arg2,  PetscInt *arg3, PetscInt n)
136: {
137:   while (n--) *arg1++ = (*arg2++ & *arg3++);
138:   return 0;
139: }

141: /***********************************ivec.c*************************************/
142: PetscInt PCTFS_ivec_sum(PetscInt *arg1,  PetscInt n)
143: {
144:   PetscInt tmp = 0;
145:   while (n--) tmp += *arg1++;
146:   return(tmp);
147: }

149: /***********************************ivec.c*************************************/
150: PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *arg1, PetscInt *arg2,  PetscInt n, ...)
151: {
152:   PetscInt i, j, type;
153:   PetscInt *arg3;
154:   va_list  ap;

156:   va_start(ap, n);
157:   arg3 = va_arg(ap, PetscInt*);
158:   va_end(ap);

160:   /* LATER: if we're really motivated we can sort and then unsort */
161:   for (i=0; i<n;) {
162:     /* clump 'em for now */
163:     j    =i+1;
164:     type = arg3[i];
165:     while ((j<n)&&(arg3[j]==type)) j++;

167:     /* how many together */
168:     j -= i;

170:     /* call appropriate ivec function */
171:     if (type == GL_MAX)        PCTFS_ivec_max(arg1,arg2,j);
172:     else if (type == GL_MIN)   PCTFS_ivec_min(arg1,arg2,j);
173:     else if (type == GL_MULT)  PCTFS_ivec_mult(arg1,arg2,j);
174:     else if (type == GL_ADD)   PCTFS_ivec_add(arg1,arg2,j);
175:     else if (type == GL_B_XOR) PCTFS_ivec_xor(arg1,arg2,j);
176:     else if (type == GL_B_OR)  PCTFS_ivec_or(arg1,arg2,j);
177:     else if (type == GL_B_AND) PCTFS_ivec_and(arg1,arg2,j);
178:     else if (type == GL_L_XOR) PCTFS_ivec_lxor(arg1,arg2,j);
179:     else if (type == GL_L_OR)  PCTFS_ivec_lor(arg1,arg2,j);
180:     else if (type == GL_L_AND) PCTFS_ivec_land(arg1,arg2,j);
181:     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_ivec_non_uniform()!");

183:     arg1+=j; arg2+=j; i+=j;
184:   }
185:   return 0;
186: }

188: /***********************************ivec.c*************************************/
189: vfp PCTFS_ivec_fct_addr(PetscInt type)
190: {
191:   if (type == NON_UNIFORM)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_non_uniform);
192:   else if (type == GL_MAX)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_max);
193:   else if (type == GL_MIN)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_min);
194:   else if (type == GL_MULT)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_mult);
195:   else if (type == GL_ADD)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_add);
196:   else if (type == GL_B_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_xor);
197:   else if (type == GL_B_OR)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_or);
198:   else if (type == GL_B_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_and);
199:   else if (type == GL_L_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lxor);
200:   else if (type == GL_L_OR)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lor);
201:   else if (type == GL_L_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_land);

203:   /* catch all ... not good if we get here */
204:   return(NULL);
205: }

207: /******************************************************************************/
208: PetscErrorCode PCTFS_ivec_sort(PetscInt *ar,  PetscInt size)
209: {
210:   PetscInt *pi, *pj, temp;
211:   PetscInt **top_a = (PetscInt**) offset_stack;
212:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;

214:   /* we're really interested in the offset of the last element */
215:   /* ==> length of the list is now size + 1                    */
216:   size--;

218:   /* do until we're done ... return when stack is exhausted */
219:   for (;;) {
220:     /* if list is large enough use quicksort partition exchange code */
221:     if (size > SORT_OPT) {
222:       /* start up pointer at element 1 and down at size     */
223:       pi = ar+1;
224:       pj = ar+size;

226:       /* find middle element in list and swap w/ element 1 */
227:       SWAP(*(ar+(size>>1)),*pi)

229:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
230:       /* note ==> pivot_value in index 0                   */
231:       if (*pi > *pj) { SWAP(*pi,*pj) }
232:       if (*ar > *pj) { SWAP(*ar,*pj) }
233:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) }

235:       /* partition about pivot_value ...                              */
236:       /* note lists of length 2 are not guaranteed to be sorted */
237:       for (;;) {
238:         /* walk up ... and down ... swap if equal to pivot! */
239:         do pi++; while (*pi<*ar);
240:         do pj--; while (*pj>*ar);

242:         /* if we've crossed we're done */
243:         if (pj<pi) break;

245:         /* else swap */
246:         SWAP(*pi,*pj)
247:       }

249:       /* place pivot_value in it's correct location */
250:       SWAP(*ar,*pj)

252:       /* test stack_size to see if we've exhausted our stack */

255:       /* push right hand child iff length > 1 */
256:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
257:         *(top_a++) = pi;
258:         size      -= *top_s+2;
259:         top_s++;
260:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
261:       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
262:         ar   = *(--top_a);
263:         size = *(--top_s);
264:       }
265:     } else { /* else sort small list directly then pop another off stack */

267:       /* insertion sort for bottom */
268:       for (pj=ar+1; pj<=ar+size; pj++) {
269:         temp = *pj;
270:         for (pi=pj-1; pi>=ar; pi--) {
271:           if (*pi <= temp) break;
272:           *(pi+1)=*pi;
273:         }
274:         *(pi+1)=temp;
275:       }

277:       /* check to see if stack is exhausted ==> DONE */
278:       if (top_s==bottom_s) return 0;

280:       /* else pop another list from the stack */
281:       ar   = *(--top_a);
282:       size = *(--top_s);
283:     }
284:   }
285: }

287: /******************************************************************************/
288: PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *ar,  PetscInt *ar2,  PetscInt size)
289: {
290:   PetscInt *pi, *pj, temp, temp2;
291:   PetscInt **top_a = (PetscInt**)offset_stack;
292:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;
293:   PetscInt *pi2, *pj2;
294:   PetscInt mid;

296:   /* we're really interested in the offset of the last element */
297:   /* ==> length of the list is now size + 1                    */
298:   size--;

300:   /* do until we're done ... return when stack is exhausted */
301:   for (;;) {

303:     /* if list is large enough use quicksort partition exchange code */
304:     if (size > SORT_OPT) {

306:       /* start up pointer at element 1 and down at size     */
307:       mid = size>>1;
308:       pi  = ar+1;
309:       pj  = ar+mid;
310:       pi2 = ar2+1;
311:       pj2 = ar2+mid;

313:       /* find middle element in list and swap w/ element 1 */
314:       SWAP(*pi,*pj)
315:       SWAP(*pi2,*pj2)

317:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
318:       /* note ==> pivot_value in index 0                   */
319:       pj  = ar+size;
320:       pj2 = ar2+size;
321:       if (*pi > *pj) { SWAP(*pi,*pj) SWAP(*pi2,*pj2) }
322:       if (*ar > *pj) { SWAP(*ar,*pj) SWAP(*ar2,*pj2) }
323:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1)) }

325:       /* partition about pivot_value ...                              */
326:       /* note lists of length 2 are not guaranteed to be sorted */
327:       for (;;) {
328:         /* walk up ... and down ... swap if equal to pivot! */
329:         do { pi++; pi2++; } while (*pi<*ar);
330:         do { pj--; pj2--; } while (*pj>*ar);

332:         /* if we've crossed we're done */
333:         if (pj<pi) break;

335:         /* else swap */
336:         SWAP(*pi,*pj)
337:         SWAP(*pi2,*pj2)
338:       }

340:       /* place pivot_value in it's correct location */
341:       SWAP(*ar,*pj)
342:       SWAP(*ar2,*pj2)

344:       /* test stack_size to see if we've exhausted our stack */

347:       /* push right hand child iff length > 1 */
348:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
349:         *(top_a++) = pi;
350:         *(top_a++) = pi2;
351:         size      -= *top_s+2;
352:         top_s++;
353:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
354:       else {  /* might as well pop - note NR_OPT >=2 ==> we're ok! */
355:         ar2  = *(--top_a);
356:         ar   = *(--top_a);
357:         size = *(--top_s);
358:       }
359:     } else { /* else sort small list directly then pop another off stack */

361:       /* insertion sort for bottom */
362:       for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
363:         temp  = *pj;
364:         temp2 = *pj2;
365:         for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
366:           if (*pi <= temp) break;
367:           *(pi+1) =*pi;
368:           *(pi2+1)=*pi2;
369:         }
370:         *(pi+1) =temp;
371:         *(pi2+1)=temp2;
372:       }

374:       /* check to see if stack is exhausted ==> DONE */
375:       if (top_s==bottom_s) return 0;

377:       /* else pop another list from the stack */
378:       ar2  = *(--top_a);
379:       ar   = *(--top_a);
380:       size = *(--top_s);
381:     }
382:   }
383: }

385: /******************************************************************************/
386: PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *ar,  PetscInt **ar2, PetscInt size)
387: {
388:   PetscInt *pi, *pj, temp, *ptr;
389:   PetscInt **top_a = (PetscInt**)offset_stack;
390:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;
391:   PetscInt **pi2, **pj2;
392:   PetscInt mid;

394:   /* we're really interested in the offset of the last element */
395:   /* ==> length of the list is now size + 1                    */
396:   size--;

398:   /* do until we're done ... return when stack is exhausted */
399:   for (;;) {

401:     /* if list is large enough use quicksort partition exchange code */
402:     if (size > SORT_OPT) {

404:       /* start up pointer at element 1 and down at size     */
405:       mid = size>>1;
406:       pi  = ar+1;
407:       pj  = ar+mid;
408:       pi2 = ar2+1;
409:       pj2 = ar2+mid;

411:       /* find middle element in list and swap w/ element 1 */
412:       SWAP(*pi,*pj)
413:       P_SWAP(*pi2,*pj2)

415:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
416:       /* note ==> pivot_value in index 0                   */
417:       pj  = ar+size;
418:       pj2 = ar2+size;
419:       if (*pi > *pj) { SWAP(*pi,*pj) P_SWAP(*pi2,*pj2) }
420:       if (*ar > *pj) { SWAP(*ar,*pj) P_SWAP(*ar2,*pj2) }
421:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1)) }

423:       /* partition about pivot_value ...                              */
424:       /* note lists of length 2 are not guaranteed to be sorted */
425:       for (;;) {

427:         /* walk up ... and down ... swap if equal to pivot! */
428:         do {pi++; pi2++;} while (*pi<*ar);
429:         do {pj--; pj2--;} while (*pj>*ar);

431:         /* if we've crossed we're done */
432:         if (pj<pi) break;

434:         /* else swap */
435:         SWAP(*pi,*pj)
436:         P_SWAP(*pi2,*pj2)
437:       }

439:       /* place pivot_value in it's correct location */
440:       SWAP(*ar,*pj)
441:       P_SWAP(*ar2,*pj2)

443:       /* test stack_size to see if we've exhausted our stack */

446:       /* push right hand child iff length > 1 */
447:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
448:         *(top_a++) = pi;
449:         *(top_a++) = (PetscInt*) pi2;
450:         size      -= *top_s+2;
451:         top_s++;
452:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
453:       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
454:         ar2  = (PetscInt**) *(--top_a);
455:         ar   = *(--top_a);
456:         size = *(--top_s);
457:       }
458:     } else  { /* else sort small list directly then pop another off stack */
459:       /* insertion sort for bottom */
460:       for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
461:         temp = *pj;
462:         ptr  = *pj2;
463:         for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
464:           if (*pi <= temp) break;
465:           *(pi+1) =*pi;
466:           *(pi2+1)=*pi2;
467:         }
468:         *(pi+1) =temp;
469:         *(pi2+1)=ptr;
470:       }

472:       /* check to see if stack is exhausted ==> DONE */
473:       if (top_s==bottom_s) return 0;

475:       /* else pop another list from the stack */
476:       ar2  = (PetscInt**)*(--top_a);
477:       ar   = *(--top_a);
478:       size = *(--top_s);
479:     }
480:   }
481: }

483: /******************************************************************************/
484: PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
485: {
486:   if (type == SORT_INTEGER) {
487:     if (ar2) PCTFS_ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);
488:     else PCTFS_ivec_sort((PetscInt*)ar1,size);
489:   } else if (type == SORT_INT_PTR) {
490:     if (ar2) PCTFS_ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt**)ar2,size);
491:     else PCTFS_ivec_sort((PetscInt*)ar1,size);
492:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_SMI_sort only does SORT_INTEGER!");
493:   return 0;
494: }

496: /***********************************ivec.c*************************************/
497: PetscInt PCTFS_ivec_linear_search(PetscInt item,  PetscInt *list,  PetscInt n)
498: {
499:   PetscInt tmp = n-1;

501:   while (n--) {
502:     if (*list++ == item) return(tmp-n);
503:   }
504:   return(-1);
505: }

507: /***********************************ivec.c*************************************/
508: PetscInt PCTFS_ivec_binary_search(PetscInt item,  PetscInt *list,  PetscInt rh)
509: {
510:   PetscInt mid, lh=0;

512:   rh--;
513:   while (lh<=rh) {
514:     mid = (lh+rh)>>1;
515:     if (*(list+mid) == item) return(mid);
516:     if (*(list+mid) > item) rh = mid-1;
517:     else lh = mid+1;
518:   }
519:   return(-1);
520: }

522: /*********************************ivec.c*************************************/
523: PetscErrorCode PCTFS_rvec_copy(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
524: {
525:   while (n--) *arg1++ = *arg2++;
526:   return 0;
527: }

529: /*********************************ivec.c*************************************/
530: PetscErrorCode PCTFS_rvec_zero(PetscScalar *arg1,  PetscInt n)
531: {
532:   while (n--) *arg1++ = 0.0;
533:   return 0;
534: }

536: /***********************************ivec.c*************************************/
537: PetscErrorCode PCTFS_rvec_one(PetscScalar *arg1,  PetscInt n)
538: {
539:   while (n--) *arg1++ = 1.0;
540:   return 0;
541: }

543: /***********************************ivec.c*************************************/
544: PetscErrorCode PCTFS_rvec_set(PetscScalar *arg1,  PetscScalar arg2,  PetscInt n)
545: {
546:   while (n--) *arg1++ = arg2;
547:   return 0;
548: }

550: /***********************************ivec.c*************************************/
551: PetscErrorCode PCTFS_rvec_scale(PetscScalar *arg1,  PetscScalar arg2,  PetscInt n)
552: {
553:   while (n--) *arg1++ *= arg2;
554:   return 0;
555: }

557: /*********************************ivec.c*************************************/
558: PetscErrorCode PCTFS_rvec_add(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
559: {
560:   while (n--) *arg1++ += *arg2++;
561:   return 0;
562: }

564: /*********************************ivec.c*************************************/
565: PetscErrorCode PCTFS_rvec_mult(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
566: {
567:   while (n--) *arg1++ *= *arg2++;
568:   return 0;
569: }

571: /*********************************ivec.c*************************************/
572: PetscErrorCode PCTFS_rvec_max(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
573: {
574:   while (n--) {
575:     *arg1 = PetscMax(*arg1,*arg2);
576:     arg1++;
577:     arg2++;
578:   }
579:   return 0;
580: }

582: /*********************************ivec.c*************************************/
583: PetscErrorCode PCTFS_rvec_max_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
584: {
585:   while (n--) {
586:     *arg1 = MAX_FABS(*arg1,*arg2);
587:     arg1++;
588:     arg2++;
589:   }
590:   return 0;
591: }

593: /*********************************ivec.c*************************************/
594: PetscErrorCode PCTFS_rvec_min(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
595: {
596:   while (n--) {
597:     *arg1 = PetscMin(*arg1,*arg2);
598:     arg1++;
599:     arg2++;
600:   }
601:   return 0;
602: }

604: /*********************************ivec.c*************************************/
605: PetscErrorCode PCTFS_rvec_min_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
606: {
607:   while (n--) {
608:     *arg1 = MIN_FABS(*arg1,*arg2);
609:     arg1++;
610:     arg2++;
611:   }
612:   return 0;
613: }

615: /*********************************ivec.c*************************************/
616: PetscErrorCode PCTFS_rvec_exists(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
617: {
618:   while (n--) {
619:     *arg1 = EXISTS(*arg1,*arg2);
620:     arg1++;
621:     arg2++;
622:   }
623:   return 0;
624: }

626: /***********************************ivec.c*************************************/
627: PetscErrorCode PCTFS_rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2,  PetscInt n,  PetscInt *arg3)
628: {
629:   PetscInt i, j, type;

631:   /* LATER: if we're really motivated we can sort and then unsort */
632:   for (i=0; i<n;) {

634:     /* clump 'em for now */
635:     j    =i+1;
636:     type = arg3[i];
637:     while ((j<n)&&(arg3[j]==type)) j++;

639:     /* how many together */
640:     j -= i;

642:     /* call appropriate ivec function */
643:     if (type == GL_MAX)          PCTFS_rvec_max(arg1,arg2,j);
644:     else if (type == GL_MIN)     PCTFS_rvec_min(arg1,arg2,j);
645:     else if (type == GL_MULT)    PCTFS_rvec_mult(arg1,arg2,j);
646:     else if (type == GL_ADD)     PCTFS_rvec_add(arg1,arg2,j);
647:     else if (type == GL_MAX_ABS) PCTFS_rvec_max_abs(arg1,arg2,j);
648:     else if (type == GL_MIN_ABS) PCTFS_rvec_min_abs(arg1,arg2,j);
649:     else if (type == GL_EXISTS)  PCTFS_rvec_exists(arg1,arg2,j);
650:     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_rvec_non_uniform()!");

652:     arg1+=j; arg2+=j; i+=j;
653:   }
654:   return 0;
655: }

657: /***********************************ivec.c*************************************/
658: vfp PCTFS_rvec_fct_addr(PetscInt type)
659: {
660:   if (type == NON_UNIFORM)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_non_uniform);
661:   else if (type == GL_MAX)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_max);
662:   else if (type == GL_MIN)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_min);
663:   else if (type == GL_MULT)    return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_mult);
664:   else if (type == GL_ADD)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_add);
665:   else if (type == GL_MAX_ABS) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_max_abs);
666:   else if (type == GL_MIN_ABS) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_min_abs);
667:   else if (type == GL_EXISTS)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_exists);

669:   /* catch all ... not good if we get here */
670:   return(NULL);
671: }