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: {
38: while (n--) *arg1++ = 0;
39: return(0);
40: }
42: /***********************************ivec.c*************************************/
43: PetscErrorCode PCTFS_ivec_set(PetscInt *arg1, PetscInt arg2, PetscInt n)
44: {
46: while (n--) *arg1++ = arg2;
47: return(0);
48: }
50: /***********************************ivec.c*************************************/
51: PetscErrorCode PCTFS_ivec_max(PetscInt *arg1, PetscInt *arg2, PetscInt n)
52: {
54: while (n--) { *arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++; }
55: return(0);
56: }
58: /***********************************ivec.c*************************************/
59: PetscErrorCode PCTFS_ivec_min(PetscInt *arg1, PetscInt *arg2, PetscInt n)
60: {
62: while (n--) {
63: *(arg1) = PetscMin(*arg1,*arg2);
64: arg1++;
65: arg2++;
66: }
67: return(0);
68: }
70: /***********************************ivec.c*************************************/
71: PetscErrorCode PCTFS_ivec_mult(PetscInt *arg1, PetscInt *arg2, PetscInt n)
72: {
74: while (n--) *arg1++ *= *arg2++;
75: return(0);
76: }
78: /***********************************ivec.c*************************************/
79: PetscErrorCode PCTFS_ivec_add(PetscInt *arg1, PetscInt *arg2, PetscInt n)
80: {
82: while (n--) *arg1++ += *arg2++;
83: return(0);
84: }
86: /***********************************ivec.c*************************************/
87: PetscErrorCode PCTFS_ivec_lxor(PetscInt *arg1, PetscInt *arg2, PetscInt n)
88: {
90: while (n--) {
91: *arg1=((*arg1 || *arg2) && !(*arg1 && *arg2));
92: arg1++;
93: arg2++;
94: }
95: return(0);
96: }
98: /***********************************ivec.c*************************************/
99: PetscErrorCode PCTFS_ivec_xor(PetscInt *arg1, PetscInt *arg2, PetscInt n)
100: {
102: while (n--) *arg1++ ^= *arg2++;
103: return(0);
104: }
106: /***********************************ivec.c*************************************/
107: PetscErrorCode PCTFS_ivec_or(PetscInt *arg1, PetscInt *arg2, PetscInt n)
108: {
110: while (n--) *arg1++ |= *arg2++;
111: return(0);
112: }
114: /***********************************ivec.c*************************************/
115: PetscErrorCode PCTFS_ivec_lor(PetscInt *arg1, PetscInt *arg2, PetscInt n)
116: {
118: while (n--) {
119: *arg1 = (*arg1 || *arg2);
120: arg1++;
121: arg2++;
122: }
123: return(0);
124: }
126: /***********************************ivec.c*************************************/
127: PetscErrorCode PCTFS_ivec_and(PetscInt *arg1, PetscInt *arg2, PetscInt n)
128: {
130: while (n--) *arg1++ &= *arg2++;
131: return(0);
132: }
134: /***********************************ivec.c*************************************/
135: PetscErrorCode PCTFS_ivec_land(PetscInt *arg1, PetscInt *arg2, PetscInt n)
136: {
138: while (n--) {
139: *arg1 = (*arg1 && *arg2);
140: arg1++;
141: arg2++;
142: }
143: return(0);
144: }
146: /***********************************ivec.c*************************************/
147: PetscErrorCode PCTFS_ivec_and3(PetscInt *arg1, PetscInt *arg2, PetscInt *arg3, PetscInt n)
148: {
150: while (n--) *arg1++ = (*arg2++ & *arg3++);
151: return(0);
152: }
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: }
162: /***********************************ivec.c*************************************/
163: PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *arg1, PetscInt *arg2, PetscInt n, ...)
164: {
165: PetscInt i, j, type;
166: PetscInt *arg3;
167: va_list ap;
170: va_start(ap, n);
171: arg3 = va_arg(ap, PetscInt*);
172: va_end(ap);
174: /* LATER: if we're really motivated we can sort and then unsort */
175: for (i=0; i<n;) {
176: /* clump 'em for now */
177: j =i+1;
178: type = arg3[i];
179: while ((j<n)&&(arg3[j]==type)) j++;
181: /* how many together */
182: j -= i;
184: /* call appropriate ivec function */
185: if (type == GL_MAX) PCTFS_ivec_max(arg1,arg2,j);
186: else if (type == GL_MIN) PCTFS_ivec_min(arg1,arg2,j);
187: else if (type == GL_MULT) PCTFS_ivec_mult(arg1,arg2,j);
188: else if (type == GL_ADD) PCTFS_ivec_add(arg1,arg2,j);
189: else if (type == GL_B_XOR) PCTFS_ivec_xor(arg1,arg2,j);
190: else if (type == GL_B_OR) PCTFS_ivec_or(arg1,arg2,j);
191: else if (type == GL_B_AND) PCTFS_ivec_and(arg1,arg2,j);
192: else if (type == GL_L_XOR) PCTFS_ivec_lxor(arg1,arg2,j);
193: else if (type == GL_L_OR) PCTFS_ivec_lor(arg1,arg2,j);
194: else if (type == GL_L_AND) PCTFS_ivec_land(arg1,arg2,j);
195: else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_ivec_non_uniform()!");
197: arg1+=j; arg2+=j; i+=j;
198: }
199: return(0);
200: }
202: /***********************************ivec.c*************************************/
203: vfp PCTFS_ivec_fct_addr(PetscInt type)
204: {
205: if (type == NON_UNIFORM) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_non_uniform);
206: else if (type == GL_MAX) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_max);
207: else if (type == GL_MIN) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_min);
208: else if (type == GL_MULT) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_mult);
209: else if (type == GL_ADD) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_add);
210: else if (type == GL_B_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_xor);
211: else if (type == GL_B_OR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_or);
212: else if (type == GL_B_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_and);
213: else if (type == GL_L_XOR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lxor);
214: else if (type == GL_L_OR) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_lor);
215: else if (type == GL_L_AND) return((PetscErrorCode (*)(void*, void*, PetscInt, ...))&PCTFS_ivec_land);
217: /* catch all ... not good if we get here */
218: return(NULL);
219: }
221: /******************************************************************************/
222: PetscErrorCode PCTFS_ivec_sort(PetscInt *ar, PetscInt size)
223: {
224: PetscInt *pi, *pj, temp;
225: PetscInt **top_a = (PetscInt**) offset_stack;
226: PetscInt *top_s = size_stack, *bottom_s = size_stack;
229: /* we're really interested in the offset of the last element */
230: /* ==> length of the list is now size + 1 */
231: size--;
233: /* do until we're done ... return when stack is exhausted */
234: for (;;) {
235: /* if list is large enough use quicksort partition exchange code */
236: if (size > SORT_OPT) {
237: /* start up pointer at element 1 and down at size */
238: pi = ar+1;
239: pj = ar+size;
241: /* find middle element in list and swap w/ element 1 */
242: SWAP(*(ar+(size>>1)),*pi)
244: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
245: /* note ==> pivot_value in index 0 */
246: if (*pi > *pj) { SWAP(*pi,*pj) }
247: if (*ar > *pj) { SWAP(*ar,*pj) }
248: else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) }
250: /* partition about pivot_value ... */
251: /* note lists of length 2 are not guaranteed to be sorted */
252: for (;;) {
253: /* walk up ... and down ... swap if equal to pivot! */
254: do pi++; while (*pi<*ar);
255: do pj--; while (*pj>*ar);
257: /* if we've crossed we're done */
258: if (pj<pi) break;
260: /* else swap */
261: SWAP(*pi,*pj)
262: }
264: /* place pivot_value in it's correct location */
265: SWAP(*ar,*pj)
267: /* test stack_size to see if we've exhausted our stack */
268: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort() :: STACK EXHAUSTED!!!");
270: /* push right hand child iff length > 1 */
271: if ((*top_s = size-((PetscInt) (pi-ar)))) {
272: *(top_a++) = pi;
273: size -= *top_s+2;
274: top_s++;
275: } else if (size -= *top_s+2) ; /* set up for next loop iff there is something to do */
276: else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
277: ar = *(--top_a);
278: size = *(--top_s);
279: }
280: } else { /* else sort small list directly then pop another off stack */
282: /* insertion sort for bottom */
283: for (pj=ar+1; pj<=ar+size; pj++) {
284: temp = *pj;
285: for (pi=pj-1; pi>=ar; pi--) {
286: if (*pi <= temp) break;
287: *(pi+1)=*pi;
288: }
289: *(pi+1)=temp;
290: }
292: /* check to see if stack is exhausted ==> DONE */
293: if (top_s==bottom_s) return(0);
295: /* else pop another list from the stack */
296: ar = *(--top_a);
297: size = *(--top_s);
298: }
299: }
300: }
302: /******************************************************************************/
303: PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *ar, PetscInt *ar2, PetscInt size)
304: {
305: PetscInt *pi, *pj, temp, temp2;
306: PetscInt **top_a = (PetscInt**)offset_stack;
307: PetscInt *top_s = size_stack, *bottom_s = size_stack;
308: PetscInt *pi2, *pj2;
309: PetscInt mid;
312: /* we're really interested in the offset of the last element */
313: /* ==> length of the list is now size + 1 */
314: size--;
316: /* do until we're done ... return when stack is exhausted */
317: for (;;) {
319: /* if list is large enough use quicksort partition exchange code */
320: if (size > SORT_OPT) {
322: /* start up pointer at element 1 and down at size */
323: mid = size>>1;
324: pi = ar+1;
325: pj = ar+mid;
326: pi2 = ar2+1;
327: pj2 = ar2+mid;
329: /* find middle element in list and swap w/ element 1 */
330: SWAP(*pi,*pj)
331: SWAP(*pi2,*pj2)
333: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
334: /* note ==> pivot_value in index 0 */
335: pj = ar+size;
336: pj2 = ar2+size;
337: if (*pi > *pj) { SWAP(*pi,*pj) SWAP(*pi2,*pj2) }
338: if (*ar > *pj) { SWAP(*ar,*pj) SWAP(*ar2,*pj2) }
339: else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1)) }
341: /* partition about pivot_value ... */
342: /* note lists of length 2 are not guaranteed to be sorted */
343: for (;;) {
344: /* walk up ... and down ... swap if equal to pivot! */
345: do { pi++; pi2++; } while (*pi<*ar);
346: do { pj--; pj2--; } while (*pj>*ar);
348: /* if we've crossed we're done */
349: if (pj<pi) break;
351: /* else swap */
352: SWAP(*pi,*pj)
353: SWAP(*pi2,*pj2)
354: }
356: /* place pivot_value in it's correct location */
357: SWAP(*ar,*pj)
358: SWAP(*ar2,*pj2)
360: /* test stack_size to see if we've exhausted our stack */
361: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion() :: STACK EXHAUSTED!!!");
363: /* push right hand child iff length > 1 */
364: if ((*top_s = size-((PetscInt) (pi-ar)))) {
365: *(top_a++) = pi;
366: *(top_a++) = pi2;
367: size -= *top_s+2;
368: top_s++;
369: } else if (size -= *top_s+2) ; /* set up for next loop iff there is something to do */
370: else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
371: ar2 = *(--top_a);
372: ar = *(--top_a);
373: size = *(--top_s);
374: }
375: } else { /* else sort small list directly then pop another off stack */
377: /* insertion sort for bottom */
378: for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
379: temp = *pj;
380: temp2 = *pj2;
381: for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
382: if (*pi <= temp) break;
383: *(pi+1) =*pi;
384: *(pi2+1)=*pi2;
385: }
386: *(pi+1) =temp;
387: *(pi2+1)=temp2;
388: }
390: /* check to see if stack is exhausted ==> DONE */
391: if (top_s==bottom_s) return(0);
393: /* else pop another list from the stack */
394: ar2 = *(--top_a);
395: ar = *(--top_a);
396: size = *(--top_s);
397: }
398: }
399: }
401: /******************************************************************************/
402: PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *ar, PetscInt **ar2, PetscInt size)
403: {
404: PetscInt *pi, *pj, temp, *ptr;
405: PetscInt **top_a = (PetscInt**)offset_stack;
406: PetscInt *top_s = size_stack, *bottom_s = size_stack;
407: PetscInt **pi2, **pj2;
408: PetscInt mid;
411: /* we're really interested in the offset of the last element */
412: /* ==> length of the list is now size + 1 */
413: size--;
415: /* do until we're done ... return when stack is exhausted */
416: for (;;) {
418: /* if list is large enough use quicksort partition exchange code */
419: if (size > SORT_OPT) {
421: /* start up pointer at element 1 and down at size */
422: mid = size>>1;
423: pi = ar+1;
424: pj = ar+mid;
425: pi2 = ar2+1;
426: pj2 = ar2+mid;
428: /* find middle element in list and swap w/ element 1 */
429: SWAP(*pi,*pj)
430: P_SWAP(*pi2,*pj2)
432: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
433: /* note ==> pivot_value in index 0 */
434: pj = ar+size;
435: pj2 = ar2+size;
436: if (*pi > *pj) { SWAP(*pi,*pj) P_SWAP(*pi2,*pj2) }
437: if (*ar > *pj) { SWAP(*ar,*pj) P_SWAP(*ar2,*pj2) }
438: else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1)) }
440: /* partition about pivot_value ... */
441: /* note lists of length 2 are not guaranteed to be sorted */
442: for (;;) {
444: /* walk up ... and down ... swap if equal to pivot! */
445: do {pi++; pi2++;} while (*pi<*ar);
446: do {pj--; pj2--;} while (*pj>*ar);
448: /* if we've crossed we're done */
449: if (pj<pi) break;
451: /* else swap */
452: SWAP(*pi,*pj)
453: P_SWAP(*pi2,*pj2)
454: }
456: /* place pivot_value in it's correct location */
457: SWAP(*ar,*pj)
458: P_SWAP(*ar2,*pj2)
460: /* test stack_size to see if we've exhausted our stack */
461: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion_hack() :: STACK EXHAUSTED!!!");
463: /* push right hand child iff length > 1 */
464: if ((*top_s = size-((PetscInt) (pi-ar)))) {
465: *(top_a++) = pi;
466: *(top_a++) = (PetscInt*) pi2;
467: size -= *top_s+2;
468: top_s++;
469: } else if (size -= *top_s+2) ; /* set up for next loop iff there is something to do */
470: else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
471: ar2 = (PetscInt**) *(--top_a);
472: ar = *(--top_a);
473: size = *(--top_s);
474: }
475: } else { /* else sort small list directly then pop another off stack */
476: /* insertion sort for bottom */
477: for (pj=ar+1, pj2=ar2+1; pj<=ar+size; pj++,pj2++) {
478: temp = *pj;
479: ptr = *pj2;
480: for (pi=pj-1,pi2=pj2-1; pi>=ar; pi--,pi2--) {
481: if (*pi <= temp) break;
482: *(pi+1) =*pi;
483: *(pi2+1)=*pi2;
484: }
485: *(pi+1) =temp;
486: *(pi2+1)=ptr;
487: }
489: /* check to see if stack is exhausted ==> DONE */
490: if (top_s==bottom_s) return(0);
492: /* else pop another list from the stack */
493: ar2 = (PetscInt**)*(--top_a);
494: ar = *(--top_a);
495: size = *(--top_s);
496: }
497: }
498: }
500: /******************************************************************************/
501: PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
502: {
504: if (type == SORT_INTEGER) {
505: if (ar2) PCTFS_ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);
506: else PCTFS_ivec_sort((PetscInt*)ar1,size);
507: } else if (type == SORT_INT_PTR) {
508: if (ar2) PCTFS_ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt**)ar2,size);
509: else PCTFS_ivec_sort((PetscInt*)ar1,size);
510: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_SMI_sort only does SORT_INTEGER!");
511: return(0);
512: }
514: /***********************************ivec.c*************************************/
515: PetscInt PCTFS_ivec_linear_search(PetscInt item, PetscInt *list, PetscInt n)
516: {
517: 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: }