Actual source code: ivec.c
1: #define PETSCKSP_DLL
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 *ivec_copy( PetscInt *arg1, PetscInt *arg2, PetscInt n)
31: {
32: while (n--) {*arg1++ = *arg2++;}
33: return(arg1);
34: }
36: /***********************************ivec.c*************************************/
37: PetscErrorCode ivec_zero( PetscInt *arg1, PetscInt n)
38: {
40: while (n--) {*arg1++ = 0;}
41: return(0);
42: }
44: /***********************************ivec.c*************************************/
45: PetscErrorCode ivec_set( PetscInt *arg1, PetscInt arg2, PetscInt n)
46: {
48: while (n--) {*arg1++ = arg2;}
49: return(0);
50: }
52: /***********************************ivec.c*************************************/
53: PetscErrorCode 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 ivec_min( PetscInt *arg1, PetscInt *arg2, PetscInt n)
62: {
64: while (n--) {*(arg1) = PetscMin(*arg1,*arg2); arg1++; arg2++;}
65: return(0);
66: }
68: /***********************************ivec.c*************************************/
69: PetscErrorCode ivec_mult( PetscInt *arg1, PetscInt *arg2, PetscInt n)
70: {
72: while (n--) {*arg1++ *= *arg2++;}
73: return(0);
74: }
76: /***********************************ivec.c*************************************/
77: PetscErrorCode ivec_add( PetscInt *arg1, PetscInt *arg2, PetscInt n)
78: {
80: while (n--) {*arg1++ += *arg2++;}
81: return(0);
82: }
84: /***********************************ivec.c*************************************/
85: PetscErrorCode ivec_lxor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
86: {
88: while (n--) {*arg1=((*arg1 || *arg2) && !(*arg1 && *arg2)); arg1++; arg2++;}
89: return(0);
90: }
92: /***********************************ivec.c*************************************/
93: PetscErrorCode ivec_xor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
94: {
96: while (n--) {*arg1++ ^= *arg2++;}
97: return(0);
98: }
100: /***********************************ivec.c*************************************/
101: PetscErrorCode ivec_or( PetscInt *arg1, PetscInt *arg2, PetscInt n)
102: {
104: while (n--) {*arg1++ |= *arg2++;}
105: return(0);
106: }
108: /***********************************ivec.c*************************************/
109: PetscErrorCode ivec_lor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
110: {
112: while (n--) {*arg1 = (*arg1 || *arg2); arg1++; arg2++;}
113: return(0);
114: }
116: /***********************************ivec.c*************************************/
117: PetscErrorCode ivec_and( PetscInt *arg1, PetscInt *arg2, PetscInt n)
118: {
120: while (n--) {*arg1++ &= *arg2++;}
121: return(0);
122: }
124: /***********************************ivec.c*************************************/
125: PetscErrorCode ivec_land( PetscInt *arg1, PetscInt *arg2, PetscInt n)
126: {
128: while (n--) {*arg1 = (*arg1 && *arg2); arg1++; arg2++;}
129: return(0);
130: }
132: /***********************************ivec.c*************************************/
133: PetscErrorCode ivec_and3( PetscInt *arg1, PetscInt *arg2, PetscInt *arg3, PetscInt n)
134: {
136: while (n--) {*arg1++ = (*arg2++ & *arg3++);}
137: return(0);
138: }
140: /***********************************ivec.c*************************************/
141: PetscInt ivec_sum( PetscInt *arg1, PetscInt n)
142: {
143: PetscInt tmp = 0;
146: while (n--) {tmp += *arg1++;}
147: return(tmp);
148: }
150: /***********************************ivec.c*************************************/
151: PetscErrorCode ivec_non_uniform(PetscInt *arg1, PetscInt *arg2, PetscInt n, PetscInt *arg3)
152: {
153: PetscInt i, j, type;
157: /* LATER: if we're really motivated we can sort and then unsort */
158: for (i=0;i<n;)
159: {
160: /* clump 'em for now */
161: j=i+1;
162: type = arg3[i];
163: while ((j<n)&&(arg3[j]==type))
164: {j++;}
165:
166: /* how many together */
167: j -= i;
169: /* call appropriate ivec function */
170: if (type == GL_MAX)
171: {ivec_max(arg1,arg2,j);}
172: else if (type == GL_MIN)
173: {ivec_min(arg1,arg2,j);}
174: else if (type == GL_MULT)
175: {ivec_mult(arg1,arg2,j);}
176: else if (type == GL_ADD)
177: {ivec_add(arg1,arg2,j);}
178: else if (type == GL_B_XOR)
179: {ivec_xor(arg1,arg2,j);}
180: else if (type == GL_B_OR)
181: {ivec_or(arg1,arg2,j);}
182: else if (type == GL_B_AND)
183: {ivec_and(arg1,arg2,j);}
184: else if (type == GL_L_XOR)
185: {ivec_lxor(arg1,arg2,j);}
186: else if (type == GL_L_OR)
187: {ivec_lor(arg1,arg2,j);}
188: else if (type == GL_L_AND)
189: {ivec_land(arg1,arg2,j);}
190: else
191: {SETERRQ(PETSC_ERR_PLIB,"unrecognized type passed to ivec_non_uniform()!");}
193: arg1+=j; arg2+=j; i+=j;
194: }
195: return(0);
196: }
198: /***********************************ivec.c*************************************/
199: vfp ivec_fct_addr( PetscInt type)
200: {
202: if (type == NON_UNIFORM)
203: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_non_uniform);}
204: else if (type == GL_MAX)
205: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_max);}
206: else if (type == GL_MIN)
207: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_min);}
208: else if (type == GL_MULT)
209: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_mult);}
210: else if (type == GL_ADD)
211: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_add);}
212: else if (type == GL_B_XOR)
213: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_xor);}
214: else if (type == GL_B_OR)
215: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_or);}
216: else if (type == GL_B_AND)
217: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_and);}
218: else if (type == GL_L_XOR)
219: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_lxor);}
220: else if (type == GL_L_OR)
221: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_lor);}
222: else if (type == GL_L_AND)
223: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&ivec_land);}
225: /* catch all ... not good if we get here */
226: return(NULL);
227: }
229: /******************************************************************************/
230: PetscErrorCode ivec_sort( PetscInt *ar, PetscInt size)
231: {
232: PetscInt *pi, *pj, temp;
233: PetscInt **top_a = (PetscInt **) offset_stack;
234: PetscInt *top_s = size_stack, *bottom_s = size_stack;
237: /* we're really interested in the offset of the last element */
238: /* ==> length of the list is now size + 1 */
239: size--;
241: /* do until we're done ... return when stack is exhausted */
242: for (;;)
243: {
244: /* if list is large enough use quicksort partition exchange code */
245: if (size > SORT_OPT)
246: {
247: /* start up pointer at element 1 and down at size */
248: pi = ar+1;
249: pj = ar+size;
251: /* find middle element in list and swap w/ element 1 */
252: SWAP(*(ar+(size>>1)),*pi)
254: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
255: /* note ==> pivot_value in index 0 */
256: if (*pi > *pj)
257: {SWAP(*pi,*pj)}
258: if (*ar > *pj)
259: {SWAP(*ar,*pj)}
260: else if (*pi > *ar)
261: {SWAP(*(ar),*(ar+1))}
263: /* partition about pivot_value ... */
264: /* note lists of length 2 are not guaranteed to be sorted */
265: for(;;)
266: {
267: /* walk up ... and down ... swap if equal to pivot! */
268: do pi++; while (*pi<*ar);
269: do pj--; while (*pj>*ar);
271: /* if we've crossed we're done */
272: if (pj<pi) break;
274: /* else swap */
275: SWAP(*pi,*pj)
276: }
278: /* place pivot_value in it's correct location */
279: SWAP(*ar,*pj)
281: /* test stack_size to see if we've exhausted our stack */
282: if (top_s-bottom_s >= SORT_STACK)
283: {SETERRQ(PETSC_ERR_PLIB,"ivec_sort() :: STACK EXHAUSTED!!!");}
285: /* push right hand child iff length > 1 */
286: if ((*top_s = size-((PetscInt) (pi-ar))))
287: {
288: *(top_a++) = pi;
289: size -= *top_s+2;
290: top_s++;
291: }
292: /* set up for next loop iff there is something to do */
293: else if (size -= *top_s+2)
294: {;}
295: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
296: else
297: {
298: ar = *(--top_a);
299: size = *(--top_s);
300: }
301: }
303: /* else sort small list directly then pop another off stack */
304: else
305: {
306: /* insertion sort for bottom */
307: for (pj=ar+1;pj<=ar+size;pj++)
308: {
309: temp = *pj;
310: for (pi=pj-1;pi>=ar;pi--)
311: {
312: if (*pi <= temp) break;
313: *(pi+1)=*pi;
314: }
315: *(pi+1)=temp;
316: }
318: /* check to see if stack is exhausted ==> DONE */
319: if (top_s==bottom_s) return(0);
320:
321: /* else pop another list from the stack */
322: ar = *(--top_a);
323: size = *(--top_s);
324: }
325: }
326: return(0);
327: }
329: /******************************************************************************/
330: PetscErrorCode ivec_sort_companion( PetscInt *ar, PetscInt *ar2, PetscInt size)
331: {
332: PetscInt *pi, *pj, temp, temp2;
333: PetscInt **top_a = (PetscInt **)offset_stack;
334: PetscInt *top_s = size_stack, *bottom_s = size_stack;
335: PetscInt *pi2, *pj2;
336: PetscInt mid;
339: /* we're really interested in the offset of the last element */
340: /* ==> length of the list is now size + 1 */
341: size--;
343: /* do until we're done ... return when stack is exhausted */
344: for (;;)
345: {
346: /* if list is large enough use quicksort partition exchange code */
347: if (size > SORT_OPT)
348: {
349: /* start up pointer at element 1 and down at size */
350: mid = size>>1;
351: pi = ar+1;
352: pj = ar+mid;
353: pi2 = ar2+1;
354: pj2 = ar2+mid;
356: /* find middle element in list and swap w/ element 1 */
357: SWAP(*pi,*pj)
358: SWAP(*pi2,*pj2)
360: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
361: /* note ==> pivot_value in index 0 */
362: pj = ar+size;
363: pj2 = ar2+size;
364: if (*pi > *pj)
365: {SWAP(*pi,*pj) SWAP(*pi2,*pj2)}
366: if (*ar > *pj)
367: {SWAP(*ar,*pj) SWAP(*ar2,*pj2)}
368: else if (*pi > *ar)
369: {SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1))}
371: /* partition about pivot_value ... */
372: /* note lists of length 2 are not guaranteed to be sorted */
373: for(;;)
374: {
375: /* walk up ... and down ... swap if equal to pivot! */
376: do {pi++; pi2++;} while (*pi<*ar);
377: do {pj--; pj2--;} while (*pj>*ar);
379: /* if we've crossed we're done */
380: if (pj<pi) break;
382: /* else swap */
383: SWAP(*pi,*pj)
384: SWAP(*pi2,*pj2)
385: }
387: /* place pivot_value in it's correct location */
388: SWAP(*ar,*pj)
389: SWAP(*ar2,*pj2)
391: /* test stack_size to see if we've exhausted our stack */
392: if (top_s-bottom_s >= SORT_STACK)
393: {SETERRQ(PETSC_ERR_PLIB,"ivec_sort_companion() :: STACK EXHAUSTED!!!");}
395: /* push right hand child iff length > 1 */
396: if ((*top_s = size-((PetscInt) (pi-ar))))
397: {
398: *(top_a++) = pi;
399: *(top_a++) = pi2;
400: size -= *top_s+2;
401: top_s++;
402: }
403: /* set up for next loop iff there is something to do */
404: else if (size -= *top_s+2)
405: {;}
406: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
407: else
408: {
409: ar2 = *(--top_a);
410: ar = *(--top_a);
411: size = *(--top_s);
412: }
413: }
415: /* else sort small list directly then pop another off stack */
416: else
417: {
418: /* insertion sort for bottom */
419: for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++)
420: {
421: temp = *pj;
422: temp2 = *pj2;
423: for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--)
424: {
425: if (*pi <= temp) break;
426: *(pi+1)=*pi;
427: *(pi2+1)=*pi2;
428: }
429: *(pi+1)=temp;
430: *(pi2+1)=temp2;
431: }
433: /* check to see if stack is exhausted ==> DONE */
434: if (top_s==bottom_s) return(0);
435:
436: /* else pop another list from the stack */
437: ar2 = *(--top_a);
438: ar = *(--top_a);
439: size = *(--top_s);
440: }
441: }
442: return(0);
443: }
445: /******************************************************************************/
446: PetscErrorCode ivec_sort_companion_hack( PetscInt *ar, PetscInt **ar2, PetscInt size)
447: {
448: PetscInt *pi, *pj, temp, *ptr;
449: PetscInt **top_a = (PetscInt **)offset_stack;
450: PetscInt *top_s = size_stack, *bottom_s = size_stack;
451: PetscInt **pi2, **pj2;
452: PetscInt mid;
455: /* we're really interested in the offset of the last element */
456: /* ==> length of the list is now size + 1 */
457: size--;
459: /* do until we're done ... return when stack is exhausted */
460: for (;;)
461: {
462: /* if list is large enough use quicksort partition exchange code */
463: if (size > SORT_OPT)
464: {
465: /* start up pointer at element 1 and down at size */
466: mid = size>>1;
467: pi = ar+1;
468: pj = ar+mid;
469: pi2 = ar2+1;
470: pj2 = ar2+mid;
472: /* find middle element in list and swap w/ element 1 */
473: SWAP(*pi,*pj)
474: P_SWAP(*pi2,*pj2)
476: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
477: /* note ==> pivot_value in index 0 */
478: pj = ar+size;
479: pj2 = ar2+size;
480: if (*pi > *pj)
481: {SWAP(*pi,*pj) P_SWAP(*pi2,*pj2)}
482: if (*ar > *pj)
483: {SWAP(*ar,*pj) P_SWAP(*ar2,*pj2)}
484: else if (*pi > *ar)
485: {SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1))}
487: /* partition about pivot_value ... */
488: /* note lists of length 2 are not guaranteed to be sorted */
489: for(;;)
490: {
491: /* walk up ... and down ... swap if equal to pivot! */
492: do {pi++; pi2++;} while (*pi<*ar);
493: do {pj--; pj2--;} while (*pj>*ar);
495: /* if we've crossed we're done */
496: if (pj<pi) break;
498: /* else swap */
499: SWAP(*pi,*pj)
500: P_SWAP(*pi2,*pj2)
501: }
503: /* place pivot_value in it's correct location */
504: SWAP(*ar,*pj)
505: P_SWAP(*ar2,*pj2)
507: /* test stack_size to see if we've exhausted our stack */
508: if (top_s-bottom_s >= SORT_STACK)
509: {SETERRQ(PETSC_ERR_PLIB,"ivec_sort_companion_hack() :: STACK EXHAUSTED!!!");}
511: /* push right hand child iff length > 1 */
512: if ((*top_s = size-((PetscInt) (pi-ar))))
513: {
514: *(top_a++) = pi;
515: *(top_a++) = (PetscInt*) pi2;
516: size -= *top_s+2;
517: top_s++;
518: }
519: /* set up for next loop iff there is something to do */
520: else if (size -= *top_s+2)
521: {;}
522: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
523: else
524: {
525: ar2 = (PetscInt **) *(--top_a);
526: ar = *(--top_a);
527: size = *(--top_s);
528: }
529: }
531: /* else sort small list directly then pop another off stack */
532: else
533: {
534: /* insertion sort for bottom */
535: for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++)
536: {
537: temp = *pj;
538: ptr = *pj2;
539: for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--)
540: {
541: if (*pi <= temp) break;
542: *(pi+1)=*pi;
543: *(pi2+1)=*pi2;
544: }
545: *(pi+1)=temp;
546: *(pi2+1)=ptr;
547: }
549: /* check to see if stack is exhausted ==> DONE */
550: if (top_s==bottom_s) return(0);
551:
552: /* else pop another list from the stack */
553: ar2 = (PetscInt **)*(--top_a);
554: ar = *(--top_a);
555: size = *(--top_s);
556: }
557: }
558: return(0);
559: }
561: /******************************************************************************/
562: PetscErrorCode SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
563: {
565: if (type == SORT_INTEGER)
566: {
567: if (ar2)
568: {ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);}
569: else
570: {ivec_sort((PetscInt*)ar1,size);}
571: }
572: else if (type == SORT_INT_PTR)
573: {
574: if (ar2)
575: {ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt **)ar2,size);}
576: else
577: {ivec_sort((PetscInt*)ar1,size);}
578: }
580: else
581: {
582: SETERRQ(PETSC_ERR_PLIB,"SMI_sort only does SORT_INTEGER!");
583: }
584: return(0);
585: }
587: /***********************************ivec.c*************************************/
588: PetscInt ivec_linear_search( PetscInt item, PetscInt *list, PetscInt n)
589: {
590: PetscInt tmp = n-1;
592: while (n--) {if (*list++ == item) {return(tmp-n);}}
593: return(-1);
594: }
596: /***********************************ivec.c*************************************/
597: PetscInt ivec_binary_search( PetscInt item, PetscInt *list, PetscInt rh)
598: {
599: PetscInt mid, lh=0;
601: rh--;
602: while (lh<=rh)
603: {
604: mid = (lh+rh)>>1;
605: if (*(list+mid) == item)
606: {return(mid);}
607: if (*(list+mid) > item)
608: {rh = mid-1;}
609: else
610: {lh = mid+1;}
611: }
612: return(-1);
613: }
615: /*********************************ivec.c*************************************/
616: PetscErrorCode rvec_copy( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
617: {
619: while (n--) {*arg1++ = *arg2++;}
620: return(0);
621: }
623: /*********************************ivec.c*************************************/
624: PetscErrorCode rvec_zero( PetscScalar *arg1, PetscInt n)
625: {
627: while (n--) {*arg1++ = 0.0;}
628: return(0);
629: }
631: /***********************************ivec.c*************************************/
632: PetscErrorCode rvec_one( PetscScalar *arg1, PetscInt n)
633: {
635: while (n--) {*arg1++ = 1.0;}
636: return(0);
637: }
639: /***********************************ivec.c*************************************/
640: PetscErrorCode rvec_set( PetscScalar *arg1, PetscScalar arg2, PetscInt n)
641: {
643: while (n--) {*arg1++ = arg2;}
644: return(0);
645: }
647: /***********************************ivec.c*************************************/
648: PetscErrorCode rvec_scale( PetscScalar *arg1, PetscScalar arg2, PetscInt n)
649: {
651: while (n--) {*arg1++ *= arg2;}
652: return(0);
653: }
655: /*********************************ivec.c*************************************/
656: PetscErrorCode rvec_add( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
657: {
659: while (n--) {*arg1++ += *arg2++;}
660: return(0);
661: }
663: /*********************************ivec.c*************************************/
664: PetscErrorCode rvec_mult( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
665: {
667: while (n--) {*arg1++ *= *arg2++;}
668: return(0);
669: }
671: /*********************************ivec.c*************************************/
672: PetscErrorCode rvec_max( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
673: {
675: while (n--) {*arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++;}
676: return(0);
677: }
679: /*********************************ivec.c*************************************/
680: PetscErrorCode rvec_max_abs( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
681: {
683: while (n--) {*arg1 = MAX_FABS(*arg1,*arg2); arg1++; arg2++;}
684: return(0);
685: }
687: /*********************************ivec.c*************************************/
688: PetscErrorCode rvec_min( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
689: {
691: while (n--) {*arg1 = PetscMin(*arg1,*arg2); arg1++; arg2++;}
692: return(0);
693: }
695: /*********************************ivec.c*************************************/
696: PetscErrorCode rvec_min_abs( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
697: {
699: while (n--) {*arg1 = MIN_FABS(*arg1,*arg2); arg1++; arg2++;}
700: return(0);
701: }
703: /*********************************ivec.c*************************************/
704: PetscErrorCode rvec_exists( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
705: {
707: while (n--) {*arg1 = EXISTS(*arg1,*arg2); arg1++; arg2++;}
708: return(0);
709: }
711: /***********************************ivec.c*************************************/
712: PetscErrorCode rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2, PetscInt n, PetscInt *arg3)
713: {
714: PetscInt i, j, type;
717: /* LATER: if we're really motivated we can sort and then unsort */
718: for (i=0;i<n;)
719: {
720: /* clump 'em for now */
721: j=i+1;
722: type = arg3[i];
723: while ((j<n)&&(arg3[j]==type))
724: {j++;}
725:
726: /* how many together */
727: j -= i;
729: /* call appropriate ivec function */
730: if (type == GL_MAX)
731: {rvec_max(arg1,arg2,j);}
732: else if (type == GL_MIN)
733: {rvec_min(arg1,arg2,j);}
734: else if (type == GL_MULT)
735: {rvec_mult(arg1,arg2,j);}
736: else if (type == GL_ADD)
737: {rvec_add(arg1,arg2,j);}
738: else if (type == GL_MAX_ABS)
739: {rvec_max_abs(arg1,arg2,j);}
740: else if (type == GL_MIN_ABS)
741: {rvec_min_abs(arg1,arg2,j);}
742: else if (type == GL_EXISTS)
743: {rvec_exists(arg1,arg2,j);}
744: else
745: {SETERRQ(PETSC_ERR_PLIB,"unrecognized type passed to rvec_non_uniform()!");}
747: arg1+=j; arg2+=j; i+=j;
748: }
749: return(0);
750: }
752: /***********************************ivec.c*************************************/
753: vfp rvec_fct_addr( PetscInt type)
754: {
755: if (type == NON_UNIFORM)
756: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_non_uniform);}
757: else if (type == GL_MAX)
758: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_max);}
759: else if (type == GL_MIN)
760: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_min);}
761: else if (type == GL_MULT)
762: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_mult);}
763: else if (type == GL_ADD)
764: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_add);}
765: else if (type == GL_MAX_ABS)
766: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_max_abs);}
767: else if (type == GL_MIN_ABS)
768: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_min_abs);}
769: else if (type == GL_EXISTS)
770: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&rvec_exists);}
772: /* catch all ... not good if we get here */
773: return(NULL);
774: }