Actual source code: ivec.c

petsc-3.7.3 2016-08-01
Report Typos and Errors

  3: /**********************************ivec.c**************************************

  5: Author: Henry M. Tufo III

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

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

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

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

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


 25: /* allocate an address and size stack for sorter(s) */
 26: static void     *offset_stack[2*SORT_STACK];
 27: static PetscInt size_stack[SORT_STACK];

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

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

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

 52: /***********************************ivec.c*************************************/
 53: PetscErrorCode PCTFS_ivec_max(PetscInt *arg1, PetscInt *arg2, PetscInt n)
 54: {
 56:   while (n--) { *arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++; }
 57:   return(0);
 58: }

 60: /***********************************ivec.c*************************************/
 61: PetscErrorCode PCTFS_ivec_min(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 62: {
 64:   while (n--) {
 65:     *(arg1) = PetscMin(*arg1,*arg2);
 66:     arg1++;
 67:     arg2++;
 68:   }
 69:   return(0);
 70: }

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

 80: /***********************************ivec.c*************************************/
 81: PetscErrorCode PCTFS_ivec_add(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 82: {
 84:   while (n--) *arg1++ += *arg2++;
 85:   return(0);
 86: }

 88: /***********************************ivec.c*************************************/
 89: PetscErrorCode PCTFS_ivec_lxor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
 90: {
 92:   while (n--) {
 93:     *arg1=((*arg1 || *arg2) && !(*arg1 && *arg2));
 94:     arg1++;
 95:     arg2++;
 96:   }
 97:   return(0);
 98: }

100: /***********************************ivec.c*************************************/
101: PetscErrorCode PCTFS_ivec_xor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
102: {
104:   while (n--) *arg1++ ^= *arg2++;
105:   return(0);
106: }

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

116: /***********************************ivec.c*************************************/
117: PetscErrorCode PCTFS_ivec_lor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
118: {
120:   while (n--) {
121:     *arg1 = (*arg1 || *arg2);
122:     arg1++;
123:     arg2++;
124:   }
125:   return(0);
126: }

128: /***********************************ivec.c*************************************/
129: PetscErrorCode PCTFS_ivec_and(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
130: {
132:   while (n--) *arg1++ &= *arg2++;
133:   return(0);
134: }

136: /***********************************ivec.c*************************************/
137: PetscErrorCode PCTFS_ivec_land(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
138: {
140:   while (n--) {
141:     *arg1 = (*arg1 && *arg2);
142:     arg1++;
143:     arg2++;
144:   }
145:   return(0);
146: }

148: /***********************************ivec.c*************************************/
149: PetscErrorCode PCTFS_ivec_and3(PetscInt *arg1,  PetscInt *arg2,  PetscInt *arg3, PetscInt n)
150: {
152:   while (n--) *arg1++ = (*arg2++ & *arg3++);
153:   return(0);
154: }

156: /***********************************ivec.c*************************************/
157: PetscInt PCTFS_ivec_sum(PetscInt *arg1,  PetscInt n)
158: {
159:   PetscInt tmp = 0;
160:   while (n--) tmp += *arg1++;
161:   return(tmp);
162: }

164: /***********************************ivec.c*************************************/
165: PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *arg1, PetscInt *arg2,  PetscInt n,  PetscInt *arg3)
166: {
167:   PetscInt i, j, type;

170:   /* LATER: if we're really motivated we can sort and then unsort */
171:   for (i=0; i<n; ) {
172:     /* clump 'em for now */
173:     j    =i+1;
174:     type = arg3[i];
175:     while ((j<n)&&(arg3[j]==type)) j++;

177:     /* how many together */
178:     j -= i;

180:     /* call appropriate ivec function */
181:     if (type == GL_MAX)        PCTFS_ivec_max(arg1,arg2,j);
182:     else if (type == GL_MIN)   PCTFS_ivec_min(arg1,arg2,j);
183:     else if (type == GL_MULT)  PCTFS_ivec_mult(arg1,arg2,j);
184:     else if (type == GL_ADD)   PCTFS_ivec_add(arg1,arg2,j);
185:     else if (type == GL_B_XOR) PCTFS_ivec_xor(arg1,arg2,j);
186:     else if (type == GL_B_OR)  PCTFS_ivec_or(arg1,arg2,j);
187:     else if (type == GL_B_AND) PCTFS_ivec_and(arg1,arg2,j);
188:     else if (type == GL_L_XOR) PCTFS_ivec_lxor(arg1,arg2,j);
189:     else if (type == GL_L_OR)  PCTFS_ivec_lor(arg1,arg2,j);
190:     else if (type == GL_L_AND) PCTFS_ivec_land(arg1,arg2,j);
191:     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_ivec_non_uniform()!");

193:     arg1+=j; arg2+=j; i+=j;
194:   }
195:   return(0);
196: }

198: /***********************************ivec.c*************************************/
199: vfp PCTFS_ivec_fct_addr(PetscInt type)
200: {
201:   if (type == NON_UNIFORM)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_non_uniform);
202:   else if (type == GL_MAX)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_max);
203:   else if (type == GL_MIN)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_min);
204:   else if (type == GL_MULT)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_mult);
205:   else if (type == GL_ADD)   return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_add);
206:   else if (type == GL_B_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_xor);
207:   else if (type == GL_B_OR)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_or);
208:   else if (type == GL_B_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_and);
209:   else if (type == GL_L_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lxor);
210:   else if (type == GL_L_OR)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lor);
211:   else if (type == GL_L_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_land);

213:   /* catch all ... not good if we get here */
214:   return(NULL);
215: }

217: /******************************************************************************/
218: PetscErrorCode PCTFS_ivec_sort(PetscInt *ar,  PetscInt size)
219: {
220:   PetscInt *pi, *pj, temp;
221:   PetscInt **top_a = (PetscInt**) offset_stack;
222:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;


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

229:   /* do until we're done ... return when stack is exhausted */
230:   for (;; ) {
231:     /* if list is large enough use quicksort partition exchange code */
232:     if (size > SORT_OPT) {
233:       /* start up pointer at element 1 and down at size     */
234:       pi = ar+1;
235:       pj = ar+size;

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

240:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
241:       /* note ==> pivot_value in index 0                   */
242:       if (*pi > *pj) { SWAP(*pi,*pj) }
243:       if (*ar > *pj) { SWAP(*ar,*pj) }
244:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) }

246:       /* partition about pivot_value ...                              */
247:       /* note lists of length 2 are not guaranteed to be sorted */
248:       for (;; ) {
249:         /* walk up ... and down ... swap if equal to pivot! */
250:         do pi++; while (*pi<*ar);
251:         do pj--; while (*pj>*ar);

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

256:         /* else swap */
257:         SWAP(*pi,*pj)
258:       }

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

263:       /* test stack_size to see if we've exhausted our stack */
264:       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort() :: STACK EXHAUSTED!!!");

266:       /* push right hand child iff length > 1 */
267:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
268:         *(top_a++) = pi;
269:         size      -= *top_s+2;
270:         top_s++;
271:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
272:       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
273:         ar   = *(--top_a);
274:         size = *(--top_s);
275:       }
276:     } else { /* else sort small list directly then pop another off stack */

278:       /* insertion sort for bottom */
279:       for (pj=ar+1; pj<=ar+size; pj++) {
280:         temp = *pj;
281:         for (pi=pj-1; pi>=ar; pi--) {
282:           if (*pi <= temp) break;
283:           *(pi+1)=*pi;
284:         }
285:         *(pi+1)=temp;
286:       }

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

291:       /* else pop another list from the stack */
292:       ar   = *(--top_a);
293:       size = *(--top_s);
294:     }
295:   }
296:   return(0);
297: }

299: /******************************************************************************/
300: PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *ar,  PetscInt *ar2,  PetscInt size)
301: {
302:   PetscInt *pi, *pj, temp, temp2;
303:   PetscInt **top_a = (PetscInt**)offset_stack;
304:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;
305:   PetscInt *pi2, *pj2;
306:   PetscInt mid;

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

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

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

319:       /* start up pointer at element 1 and down at size     */
320:       mid = size>>1;
321:       pi  = ar+1;
322:       pj  = ar+mid;
323:       pi2 = ar2+1;
324:       pj2 = ar2+mid;

326:       /* find middle element in list and swap w/ element 1 */
327:       SWAP(*pi,*pj)
328:       SWAP(*pi2,*pj2)

330:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
331:       /* note ==> pivot_value in index 0                   */
332:       pj  = ar+size;
333:       pj2 = ar2+size;
334:       if (*pi > *pj) { SWAP(*pi,*pj) SWAP(*pi2,*pj2) }
335:       if (*ar > *pj) { SWAP(*ar,*pj) SWAP(*ar2,*pj2) }
336:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1)) }

338:       /* partition about pivot_value ...                              */
339:       /* note lists of length 2 are not guaranteed to be sorted */
340:       for (;; ) {
341:         /* walk up ... and down ... swap if equal to pivot! */
342:         do { pi++; pi2++; } while (*pi<*ar);
343:         do { pj--; pj2--; } while (*pj>*ar);

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

348:         /* else swap */
349:         SWAP(*pi,*pj)
350:         SWAP(*pi2,*pj2)
351:       }

353:       /* place pivot_value in it's correct location */
354:       SWAP(*ar,*pj)
355:       SWAP(*ar2,*pj2)

357:       /* test stack_size to see if we've exhausted our stack */
358:       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion() :: STACK EXHAUSTED!!!");

360:       /* push right hand child iff length > 1 */
361:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
362:         *(top_a++) = pi;
363:         *(top_a++) = pi2;
364:         size      -= *top_s+2;
365:         top_s++;
366:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
367:       else {  /* might as well pop - note NR_OPT >=2 ==> we're ok! */
368:         ar2  = *(--top_a);
369:         ar   = *(--top_a);
370:         size = *(--top_s);
371:       }
372:     } else { /* else sort small list directly then pop another off stack */

374:       /* insertion sort for bottom */
375:       for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
376:         temp  = *pj;
377:         temp2 = *pj2;
378:         for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
379:           if (*pi <= temp) break;
380:           *(pi+1) =*pi;
381:           *(pi2+1)=*pi2;
382:         }
383:         *(pi+1) =temp;
384:         *(pi2+1)=temp2;
385:       }

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

390:       /* else pop another list from the stack */
391:       ar2  = *(--top_a);
392:       ar   = *(--top_a);
393:       size = *(--top_s);
394:     }
395:   }
396:   return(0);
397: }

399: /******************************************************************************/
400: PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *ar,  PetscInt **ar2, PetscInt size)
401: {
402:   PetscInt *pi, *pj, temp, *ptr;
403:   PetscInt **top_a = (PetscInt**)offset_stack;
404:   PetscInt *top_s  = size_stack, *bottom_s = size_stack;
405:   PetscInt **pi2, **pj2;
406:   PetscInt mid;

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

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

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

419:       /* start up pointer at element 1 and down at size     */
420:       mid = size>>1;
421:       pi  = ar+1;
422:       pj  = ar+mid;
423:       pi2 = ar2+1;
424:       pj2 = ar2+mid;

426:       /* find middle element in list and swap w/ element 1 */
427:       SWAP(*pi,*pj)
428:       P_SWAP(*pi2,*pj2)

430:       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
431:       /* note ==> pivot_value in index 0                   */
432:       pj  = ar+size;
433:       pj2 = ar2+size;
434:       if (*pi > *pj) { SWAP(*pi,*pj) P_SWAP(*pi2,*pj2) }
435:       if (*ar > *pj) { SWAP(*ar,*pj) P_SWAP(*ar2,*pj2) }
436:       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1)) }

438:       /* partition about pivot_value ...                              */
439:       /* note lists of length 2 are not guaranteed to be sorted */
440:       for (;; ) {

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

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

449:         /* else swap */
450:         SWAP(*pi,*pj)
451:         P_SWAP(*pi2,*pj2)
452:       }

454:       /* place pivot_value in it's correct location */
455:       SWAP(*ar,*pj)
456:       P_SWAP(*ar2,*pj2)

458:       /* test stack_size to see if we've exhausted our stack */
459:       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion_hack() :: STACK EXHAUSTED!!!");

461:       /* push right hand child iff length > 1 */
462:       if ((*top_s = size-((PetscInt) (pi-ar)))) {
463:         *(top_a++) = pi;
464:         *(top_a++) = (PetscInt*) pi2;
465:         size      -= *top_s+2;
466:         top_s++;
467:       } else if (size -= *top_s+2) ;   /* set up for next loop iff there is something to do */
468:       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
469:         ar2  = (PetscInt**) *(--top_a);
470:         ar   = *(--top_a);
471:         size = *(--top_s);
472:       }
473:     } else  { /* else sort small list directly then pop another off stack */
474:       /* insertion sort for bottom */
475:       for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
476:         temp = *pj;
477:         ptr  = *pj2;
478:         for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
479:           if (*pi <= temp) break;
480:           *(pi+1) =*pi;
481:           *(pi2+1)=*pi2;
482:         }
483:         *(pi+1) =temp;
484:         *(pi2+1)=ptr;
485:       }

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

490:       /* else pop another list from the stack */
491:       ar2  = (PetscInt**)*(--top_a);
492:       ar   = *(--top_a);
493:       size = *(--top_s);
494:     }
495:   }
496:   return(0);
497: }

499: /******************************************************************************/
500: PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
501: {
503:   if (type == SORT_INTEGER) {
504:     if (ar2) PCTFS_ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);
505:     else PCTFS_ivec_sort((PetscInt*)ar1,size);
506:   } else if (type == SORT_INT_PTR) {
507:     if (ar2) PCTFS_ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt**)ar2,size);
508:     else PCTFS_ivec_sort((PetscInt*)ar1,size);
509:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_SMI_sort only does SORT_INTEGER!");
510:   return(0);
511: }

513: /***********************************ivec.c*************************************/
514: PetscInt PCTFS_ivec_linear_search(PetscInt item,  PetscInt *list,  PetscInt n)
515: {
516:   PetscInt tmp = n-1;

519:   while (n--) {
520:     if (*list++ == item) return(tmp-n);
521:   }
522:   return(-1);
523: }

525: /***********************************ivec.c*************************************/
526: PetscInt PCTFS_ivec_binary_search(PetscInt item,  PetscInt *list,  PetscInt rh)
527: {
528:   PetscInt mid, lh=0;

530:   rh--;
531:   while (lh<=rh) {
532:     mid = (lh+rh)>>1;
533:     if (*(list+mid) == item) return(mid);
534:     if (*(list+mid) > item) rh = mid-1;
535:     else lh = mid+1;
536:   }
537:   return(-1);
538: }

540: /*********************************ivec.c*************************************/
541: PetscErrorCode PCTFS_rvec_copy(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
542: {
544:   while (n--) *arg1++ = *arg2++;
545:   return(0);
546: }

548: /*********************************ivec.c*************************************/
549: PetscErrorCode PCTFS_rvec_zero(PetscScalar *arg1,  PetscInt n)
550: {
552:   while (n--) *arg1++ = 0.0;
553:   return(0);
554: }

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

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

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

580: /*********************************ivec.c*************************************/
581: PetscErrorCode PCTFS_rvec_add(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
582: {
584:   while (n--) *arg1++ += *arg2++;
585:   return(0);
586: }

588: /*********************************ivec.c*************************************/
589: PetscErrorCode PCTFS_rvec_mult(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
590: {
592:   while (n--) *arg1++ *= *arg2++;
593:   return(0);
594: }

596: /*********************************ivec.c*************************************/
597: PetscErrorCode PCTFS_rvec_max(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
598: {
600:   while (n--) {
601:     *arg1 = PetscMax(*arg1,*arg2);
602:     arg1++;
603:     arg2++;
604:   }
605:   return(0);
606: }

608: /*********************************ivec.c*************************************/
609: PetscErrorCode PCTFS_rvec_max_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
610: {
612:   while (n--) {
613:     *arg1 = MAX_FABS(*arg1,*arg2);
614:     arg1++;
615:     arg2++;
616:   }
617:   return(0);
618: }

620: /*********************************ivec.c*************************************/
621: PetscErrorCode PCTFS_rvec_min(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
622: {
624:   while (n--) {
625:     *arg1 = PetscMin(*arg1,*arg2);
626:     arg1++;
627:     arg2++;
628:   }
629:   return(0);
630: }

632: /*********************************ivec.c*************************************/
633: PetscErrorCode PCTFS_rvec_min_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
634: {
636:   while (n--) {
637:     *arg1 = MIN_FABS(*arg1,*arg2);
638:     arg1++;
639:     arg2++;
640:   }
641:   return(0);
642: }

644: /*********************************ivec.c*************************************/
645: PetscErrorCode PCTFS_rvec_exists(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
646: {
648:   while (n--) {
649:     *arg1 = EXISTS(*arg1,*arg2);
650:     arg1++;
651:     arg2++;
652:   }
653:   return(0);
654: }

656: /***********************************ivec.c*************************************/
657: PetscErrorCode PCTFS_rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2,  PetscInt n,  PetscInt *arg3)
658: {
659:   PetscInt i, j, type;

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

665:     /* clump 'em for now */
666:     j    =i+1;
667:     type = arg3[i];
668:     while ((j<n)&&(arg3[j]==type)) j++;

670:     /* how many together */
671:     j -= i;

673:     /* call appropriate ivec function */
674:     if (type == GL_MAX)          PCTFS_rvec_max(arg1,arg2,j);
675:     else if (type == GL_MIN)     PCTFS_rvec_min(arg1,arg2,j);
676:     else if (type == GL_MULT)    PCTFS_rvec_mult(arg1,arg2,j);
677:     else if (type == GL_ADD)     PCTFS_rvec_add(arg1,arg2,j);
678:     else if (type == GL_MAX_ABS) PCTFS_rvec_max_abs(arg1,arg2,j);
679:     else if (type == GL_MIN_ABS) PCTFS_rvec_min_abs(arg1,arg2,j);
680:     else if (type == GL_EXISTS)  PCTFS_rvec_exists(arg1,arg2,j);
681:     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_rvec_non_uniform()!");

683:     arg1+=j; arg2+=j; i+=j;
684:   }
685:   return(0);
686: }

688: /***********************************ivec.c*************************************/
689: vfp PCTFS_rvec_fct_addr(PetscInt type)
690: {
691:   if (type == NON_UNIFORM)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_non_uniform);
692:   else if (type == GL_MAX)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_max);
693:   else if (type == GL_MIN)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_min);
694:   else if (type == GL_MULT)    return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_mult);
695:   else if (type == GL_ADD)     return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_add);
696:   else if (type == GL_MAX_ABS) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_max_abs);
697:   else if (type == GL_MIN_ABS) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_min_abs);
698:   else if (type == GL_EXISTS)  return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_rvec_exists);

700:   /* catch all ... not good if we get here */
701:   return(NULL);
702: }