Actual source code: ivec.c
petsc-3.7.3 2016-08-01
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: }