Actual source code: ivec.c
petsc-3.3-p7 2013-05-11
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
24: /* allocate an address and size stack for sorter(s) */
25: static void *offset_stack[2*SORT_STACK];
26: static PetscInt size_stack[SORT_STACK];
28: /***********************************ivec.c*************************************/
29: PetscInt *PCTFS_ivec_copy( PetscInt *arg1, PetscInt *arg2, PetscInt n)
30: {
31: while (n--) {*arg1++ = *arg2++;}
32: return(arg1);
33: }
35: /***********************************ivec.c*************************************/
36: PetscErrorCode PCTFS_ivec_zero( PetscInt *arg1, PetscInt n)
37: {
39: while (n--) {*arg1++ = 0;}
40: return(0);
41: }
43: /***********************************ivec.c*************************************/
44: PetscErrorCode PCTFS_ivec_set( PetscInt *arg1, PetscInt arg2, PetscInt n)
45: {
47: while (n--) {*arg1++ = arg2;}
48: return(0);
49: }
51: /***********************************ivec.c*************************************/
52: PetscErrorCode PCTFS_ivec_max( PetscInt *arg1, PetscInt *arg2, PetscInt n)
53: {
55: while (n--) {*arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++;}
56: return(0);
57: }
59: /***********************************ivec.c*************************************/
60: PetscErrorCode PCTFS_ivec_min( PetscInt *arg1, PetscInt *arg2, PetscInt n)
61: {
63: while (n--) {*(arg1) = PetscMin(*arg1,*arg2); arg1++; arg2++;}
64: return(0);
65: }
67: /***********************************ivec.c*************************************/
68: PetscErrorCode PCTFS_ivec_mult( PetscInt *arg1, PetscInt *arg2, PetscInt n)
69: {
71: while (n--) {*arg1++ *= *arg2++;}
72: return(0);
73: }
75: /***********************************ivec.c*************************************/
76: PetscErrorCode PCTFS_ivec_add( PetscInt *arg1, PetscInt *arg2, PetscInt n)
77: {
79: while (n--) {*arg1++ += *arg2++;}
80: return(0);
81: }
83: /***********************************ivec.c*************************************/
84: PetscErrorCode PCTFS_ivec_lxor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
85: {
87: while (n--) {*arg1=((*arg1 || *arg2) && !(*arg1 && *arg2)); arg1++; arg2++;}
88: return(0);
89: }
91: /***********************************ivec.c*************************************/
92: PetscErrorCode PCTFS_ivec_xor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
93: {
95: while (n--) {*arg1++ ^= *arg2++;}
96: return(0);
97: }
99: /***********************************ivec.c*************************************/
100: PetscErrorCode PCTFS_ivec_or( PetscInt *arg1, PetscInt *arg2, PetscInt n)
101: {
103: while (n--) {*arg1++ |= *arg2++;}
104: return(0);
105: }
107: /***********************************ivec.c*************************************/
108: PetscErrorCode PCTFS_ivec_lor( PetscInt *arg1, PetscInt *arg2, PetscInt n)
109: {
111: while (n--) {*arg1 = (*arg1 || *arg2); arg1++; arg2++;}
112: return(0);
113: }
115: /***********************************ivec.c*************************************/
116: PetscErrorCode PCTFS_ivec_and( PetscInt *arg1, PetscInt *arg2, PetscInt n)
117: {
119: while (n--) {*arg1++ &= *arg2++;}
120: return(0);
121: }
123: /***********************************ivec.c*************************************/
124: PetscErrorCode PCTFS_ivec_land( PetscInt *arg1, PetscInt *arg2, PetscInt n)
125: {
127: while (n--) {*arg1 = (*arg1 && *arg2); arg1++; arg2++;}
128: return(0);
129: }
131: /***********************************ivec.c*************************************/
132: PetscErrorCode PCTFS_ivec_and3( PetscInt *arg1, PetscInt *arg2, PetscInt *arg3, PetscInt n)
133: {
135: while (n--) {*arg1++ = (*arg2++ & *arg3++);}
136: return(0);
137: }
139: /***********************************ivec.c*************************************/
140: PetscInt PCTFS_ivec_sum( PetscInt *arg1, PetscInt n)
141: {
142: 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, PetscInt *arg3)
151: {
152: PetscInt i, j, type;
156: /* LATER: if we're really motivated we can sort and then unsort */
157: for (i=0;i<n;)
158: {
159: /* clump 'em for now */
160: j=i+1;
161: type = arg3[i];
162: while ((j<n)&&(arg3[j]==type))
163: {j++;}
164:
165: /* how many together */
166: j -= i;
168: /* call appropriate ivec function */
169: if (type == GL_MAX)
170: {PCTFS_ivec_max(arg1,arg2,j);}
171: else if (type == GL_MIN)
172: {PCTFS_ivec_min(arg1,arg2,j);}
173: else if (type == GL_MULT)
174: {PCTFS_ivec_mult(arg1,arg2,j);}
175: else if (type == GL_ADD)
176: {PCTFS_ivec_add(arg1,arg2,j);}
177: else if (type == GL_B_XOR)
178: {PCTFS_ivec_xor(arg1,arg2,j);}
179: else if (type == GL_B_OR)
180: {PCTFS_ivec_or(arg1,arg2,j);}
181: else if (type == GL_B_AND)
182: {PCTFS_ivec_and(arg1,arg2,j);}
183: else if (type == GL_L_XOR)
184: {PCTFS_ivec_lxor(arg1,arg2,j);}
185: else if (type == GL_L_OR)
186: {PCTFS_ivec_lor(arg1,arg2,j);}
187: else if (type == GL_L_AND)
188: {PCTFS_ivec_land(arg1,arg2,j);}
189: else
190: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_ivec_non_uniform()!");
192: arg1+=j; arg2+=j; i+=j;
193: }
194: return(0);
195: }
197: /***********************************ivec.c*************************************/
198: vfp PCTFS_ivec_fct_addr( PetscInt type)
199: {
201: if (type == NON_UNIFORM)
202: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_non_uniform);}
203: else if (type == GL_MAX)
204: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_max);}
205: else if (type == GL_MIN)
206: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_min);}
207: else if (type == GL_MULT)
208: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_mult);}
209: else if (type == GL_ADD)
210: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_add);}
211: else if (type == GL_B_XOR)
212: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_xor);}
213: else if (type == GL_B_OR)
214: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_or);}
215: else if (type == GL_B_AND)
216: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_and);}
217: else if (type == GL_L_XOR)
218: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_lxor);}
219: else if (type == GL_L_OR)
220: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_lor);}
221: else if (type == GL_L_AND)
222: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_land);}
224: /* catch all ... not good if we get here */
225: return(NULL);
226: }
228: /******************************************************************************/
229: PetscErrorCode PCTFS_ivec_sort( PetscInt *ar, PetscInt size)
230: {
231: PetscInt *pi, *pj, temp;
232: PetscInt **top_a = (PetscInt **) offset_stack;
233: PetscInt *top_s = size_stack, *bottom_s = size_stack;
236: /* we're really interested in the offset of the last element */
237: /* ==> length of the list is now size + 1 */
238: size--;
240: /* do until we're done ... return when stack is exhausted */
241: for (;;)
242: {
243: /* if list is large enough use quicksort partition exchange code */
244: if (size > SORT_OPT)
245: {
246: /* start up pointer at element 1 and down at size */
247: pi = ar+1;
248: pj = ar+size;
250: /* find middle element in list and swap w/ element 1 */
251: SWAP(*(ar+(size>>1)),*pi)
253: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
254: /* note ==> pivot_value in index 0 */
255: if (*pi > *pj)
256: {SWAP(*pi,*pj)}
257: if (*ar > *pj)
258: {SWAP(*ar,*pj)}
259: else if (*pi > *ar)
260: {SWAP(*(ar),*(ar+1))}
262: /* partition about pivot_value ... */
263: /* note lists of length 2 are not guaranteed to be sorted */
264: for(;;)
265: {
266: /* walk up ... and down ... swap if equal to pivot! */
267: do pi++; while (*pi<*ar);
268: do pj--; while (*pj>*ar);
270: /* if we've crossed we're done */
271: if (pj<pi) break;
273: /* else swap */
274: SWAP(*pi,*pj)
275: }
277: /* place pivot_value in it's correct location */
278: SWAP(*ar,*pj)
280: /* test stack_size to see if we've exhausted our stack */
281: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort() :: STACK EXHAUSTED!!!");
283: /* push right hand child iff length > 1 */
284: if ((*top_s = size-((PetscInt) (pi-ar))))
285: {
286: *(top_a++) = pi;
287: size -= *top_s+2;
288: top_s++;
289: }
290: /* set up for next loop iff there is something to do */
291: else if (size -= *top_s+2)
292: {;}
293: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
294: else
295: {
296: ar = *(--top_a);
297: size = *(--top_s);
298: }
299: }
301: /* else sort small list directly then pop another off stack */
302: else
303: {
304: /* insertion sort for bottom */
305: for (pj=ar+1;pj<=ar+size;pj++)
306: {
307: temp = *pj;
308: for (pi=pj-1;pi>=ar;pi--)
309: {
310: if (*pi <= temp) break;
311: *(pi+1)=*pi;
312: }
313: *(pi+1)=temp;
314: }
316: /* check to see if stack is exhausted ==> DONE */
317: if (top_s==bottom_s) return(0);
318:
319: /* else pop another list from the stack */
320: ar = *(--top_a);
321: size = *(--top_s);
322: }
323: }
324: return(0);
325: }
327: /******************************************************************************/
328: PetscErrorCode PCTFS_ivec_sort_companion( PetscInt *ar, PetscInt *ar2, PetscInt size)
329: {
330: PetscInt *pi, *pj, temp, temp2;
331: PetscInt **top_a = (PetscInt **)offset_stack;
332: PetscInt *top_s = size_stack, *bottom_s = size_stack;
333: PetscInt *pi2, *pj2;
334: PetscInt mid;
337: /* we're really interested in the offset of the last element */
338: /* ==> length of the list is now size + 1 */
339: size--;
341: /* do until we're done ... return when stack is exhausted */
342: for (;;)
343: {
344: /* if list is large enough use quicksort partition exchange code */
345: if (size > SORT_OPT)
346: {
347: /* start up pointer at element 1 and down at size */
348: mid = size>>1;
349: pi = ar+1;
350: pj = ar+mid;
351: pi2 = ar2+1;
352: pj2 = ar2+mid;
354: /* find middle element in list and swap w/ element 1 */
355: SWAP(*pi,*pj)
356: SWAP(*pi2,*pj2)
358: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
359: /* note ==> pivot_value in index 0 */
360: pj = ar+size;
361: pj2 = ar2+size;
362: if (*pi > *pj)
363: {SWAP(*pi,*pj) SWAP(*pi2,*pj2)}
364: if (*ar > *pj)
365: {SWAP(*ar,*pj) SWAP(*ar2,*pj2)}
366: else if (*pi > *ar)
367: {SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1))}
369: /* partition about pivot_value ... */
370: /* note lists of length 2 are not guaranteed to be sorted */
371: for(;;)
372: {
373: /* walk up ... and down ... swap if equal to pivot! */
374: do {pi++; pi2++;} while (*pi<*ar);
375: do {pj--; pj2--;} while (*pj>*ar);
377: /* if we've crossed we're done */
378: if (pj<pi) break;
380: /* else swap */
381: SWAP(*pi,*pj)
382: SWAP(*pi2,*pj2)
383: }
385: /* place pivot_value in it's correct location */
386: SWAP(*ar,*pj)
387: SWAP(*ar2,*pj2)
389: /* test stack_size to see if we've exhausted our stack */
390: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion() :: STACK EXHAUSTED!!!");
392: /* push right hand child iff length > 1 */
393: if ((*top_s = size-((PetscInt) (pi-ar))))
394: {
395: *(top_a++) = pi;
396: *(top_a++) = pi2;
397: size -= *top_s+2;
398: top_s++;
399: }
400: /* set up for next loop iff there is something to do */
401: else if (size -= *top_s+2)
402: {;}
403: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
404: else
405: {
406: ar2 = *(--top_a);
407: ar = *(--top_a);
408: size = *(--top_s);
409: }
410: }
412: /* else sort small list directly then pop another off stack */
413: else
414: {
415: /* insertion sort for bottom */
416: for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++)
417: {
418: temp = *pj;
419: temp2 = *pj2;
420: for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--)
421: {
422: if (*pi <= temp) break;
423: *(pi+1)=*pi;
424: *(pi2+1)=*pi2;
425: }
426: *(pi+1)=temp;
427: *(pi2+1)=temp2;
428: }
430: /* check to see if stack is exhausted ==> DONE */
431: if (top_s==bottom_s) return(0);
432:
433: /* else pop another list from the stack */
434: ar2 = *(--top_a);
435: ar = *(--top_a);
436: size = *(--top_s);
437: }
438: }
439: return(0);
440: }
442: /******************************************************************************/
443: PetscErrorCode PCTFS_ivec_sort_companion_hack( PetscInt *ar, PetscInt **ar2, PetscInt size)
444: {
445: PetscInt *pi, *pj, temp, *ptr;
446: PetscInt **top_a = (PetscInt **)offset_stack;
447: PetscInt *top_s = size_stack, *bottom_s = size_stack;
448: PetscInt **pi2, **pj2;
449: PetscInt mid;
452: /* we're really interested in the offset of the last element */
453: /* ==> length of the list is now size + 1 */
454: size--;
456: /* do until we're done ... return when stack is exhausted */
457: for (;;)
458: {
459: /* if list is large enough use quicksort partition exchange code */
460: if (size > SORT_OPT)
461: {
462: /* start up pointer at element 1 and down at size */
463: mid = size>>1;
464: pi = ar+1;
465: pj = ar+mid;
466: pi2 = ar2+1;
467: pj2 = ar2+mid;
469: /* find middle element in list and swap w/ element 1 */
470: SWAP(*pi,*pj)
471: P_SWAP(*pi2,*pj2)
473: /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
474: /* note ==> pivot_value in index 0 */
475: pj = ar+size;
476: pj2 = ar2+size;
477: if (*pi > *pj)
478: {SWAP(*pi,*pj) P_SWAP(*pi2,*pj2)}
479: if (*ar > *pj)
480: {SWAP(*ar,*pj) P_SWAP(*ar2,*pj2)}
481: else if (*pi > *ar)
482: {SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1))}
484: /* partition about pivot_value ... */
485: /* note lists of length 2 are not guaranteed to be sorted */
486: for(;;)
487: {
488: /* walk up ... and down ... swap if equal to pivot! */
489: do {pi++; pi2++;} while (*pi<*ar);
490: do {pj--; pj2--;} while (*pj>*ar);
492: /* if we've crossed we're done */
493: if (pj<pi) break;
495: /* else swap */
496: SWAP(*pi,*pj)
497: P_SWAP(*pi2,*pj2)
498: }
500: /* place pivot_value in it's correct location */
501: SWAP(*ar,*pj)
502: P_SWAP(*ar2,*pj2)
504: /* test stack_size to see if we've exhausted our stack */
505: if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion_hack() :: STACK EXHAUSTED!!!");
507: /* push right hand child iff length > 1 */
508: if ((*top_s = size-((PetscInt) (pi-ar))))
509: {
510: *(top_a++) = pi;
511: *(top_a++) = (PetscInt*) pi2;
512: size -= *top_s+2;
513: top_s++;
514: }
515: /* set up for next loop iff there is something to do */
516: else if (size -= *top_s+2)
517: {;}
518: /* might as well pop - note NR_OPT >=2 ==> we're ok! */
519: else
520: {
521: ar2 = (PetscInt **) *(--top_a);
522: ar = *(--top_a);
523: size = *(--top_s);
524: }
525: }
527: /* else sort small list directly then pop another off stack */
528: else
529: {
530: /* insertion sort for bottom */
531: for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++)
532: {
533: temp = *pj;
534: ptr = *pj2;
535: for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--)
536: {
537: if (*pi <= temp) break;
538: *(pi+1)=*pi;
539: *(pi2+1)=*pi2;
540: }
541: *(pi+1)=temp;
542: *(pi2+1)=ptr;
543: }
545: /* check to see if stack is exhausted ==> DONE */
546: if (top_s==bottom_s) return(0);
547:
548: /* else pop another list from the stack */
549: ar2 = (PetscInt **)*(--top_a);
550: ar = *(--top_a);
551: size = *(--top_s);
552: }
553: }
554: return(0);
555: }
557: /******************************************************************************/
558: PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
559: {
561: if (type == SORT_INTEGER) {
562: if (ar2) PCTFS_ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);
563: else PCTFS_ivec_sort((PetscInt*)ar1,size);
564: } else if (type == SORT_INT_PTR) {
565: if (ar2) PCTFS_ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt **)ar2,size);
566: else PCTFS_ivec_sort((PetscInt*)ar1,size);
567: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_SMI_sort only does SORT_INTEGER!");
568: return(0);
569: }
571: /***********************************ivec.c*************************************/
572: PetscInt PCTFS_ivec_linear_search( PetscInt item, PetscInt *list, PetscInt n)
573: {
574: PetscInt tmp = n-1;
576: while (n--) {if (*list++ == item) {return(tmp-n);}}
577: return(-1);
578: }
580: /***********************************ivec.c*************************************/
581: PetscInt PCTFS_ivec_binary_search( PetscInt item, PetscInt *list, PetscInt rh)
582: {
583: PetscInt mid, lh=0;
585: rh--;
586: while (lh<=rh)
587: {
588: mid = (lh+rh)>>1;
589: if (*(list+mid) == item)
590: {return(mid);}
591: if (*(list+mid) > item)
592: {rh = mid-1;}
593: else
594: {lh = mid+1;}
595: }
596: return(-1);
597: }
599: /*********************************ivec.c*************************************/
600: PetscErrorCode PCTFS_rvec_copy( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
601: {
603: while (n--) {*arg1++ = *arg2++;}
604: return(0);
605: }
607: /*********************************ivec.c*************************************/
608: PetscErrorCode PCTFS_rvec_zero( PetscScalar *arg1, PetscInt n)
609: {
611: while (n--) {*arg1++ = 0.0;}
612: return(0);
613: }
615: /***********************************ivec.c*************************************/
616: PetscErrorCode PCTFS_rvec_one( PetscScalar *arg1, PetscInt n)
617: {
619: while (n--) {*arg1++ = 1.0;}
620: return(0);
621: }
623: /***********************************ivec.c*************************************/
624: PetscErrorCode PCTFS_rvec_set( PetscScalar *arg1, PetscScalar arg2, PetscInt n)
625: {
627: while (n--) {*arg1++ = arg2;}
628: return(0);
629: }
631: /***********************************ivec.c*************************************/
632: PetscErrorCode PCTFS_rvec_scale( PetscScalar *arg1, PetscScalar arg2, PetscInt n)
633: {
635: while (n--) {*arg1++ *= arg2;}
636: return(0);
637: }
639: /*********************************ivec.c*************************************/
640: PetscErrorCode PCTFS_rvec_add( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
641: {
643: while (n--) {*arg1++ += *arg2++;}
644: return(0);
645: }
647: /*********************************ivec.c*************************************/
648: PetscErrorCode PCTFS_rvec_mult( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
649: {
651: while (n--) {*arg1++ *= *arg2++;}
652: return(0);
653: }
655: /*********************************ivec.c*************************************/
656: PetscErrorCode PCTFS_rvec_max( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
657: {
659: while (n--) {*arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++;}
660: return(0);
661: }
663: /*********************************ivec.c*************************************/
664: PetscErrorCode PCTFS_rvec_max_abs( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
665: {
667: while (n--) {*arg1 = MAX_FABS(*arg1,*arg2); arg1++; arg2++;}
668: return(0);
669: }
671: /*********************************ivec.c*************************************/
672: PetscErrorCode PCTFS_rvec_min( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
673: {
675: while (n--) {*arg1 = PetscMin(*arg1,*arg2); arg1++; arg2++;}
676: return(0);
677: }
679: /*********************************ivec.c*************************************/
680: PetscErrorCode PCTFS_rvec_min_abs( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
681: {
683: while (n--) {*arg1 = MIN_FABS(*arg1,*arg2); arg1++; arg2++;}
684: return(0);
685: }
687: /*********************************ivec.c*************************************/
688: PetscErrorCode PCTFS_rvec_exists( PetscScalar *arg1, PetscScalar *arg2, PetscInt n)
689: {
691: while (n--) {*arg1 = EXISTS(*arg1,*arg2); arg1++; arg2++;}
692: return(0);
693: }
695: /***********************************ivec.c*************************************/
696: PetscErrorCode PCTFS_rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2, PetscInt n, PetscInt *arg3)
697: {
698: PetscInt i, j, type;
701: /* LATER: if we're really motivated we can sort and then unsort */
702: for (i=0;i<n;)
703: {
704: /* clump 'em for now */
705: j=i+1;
706: type = arg3[i];
707: while ((j<n)&&(arg3[j]==type))
708: {j++;}
709:
710: /* how many together */
711: j -= i;
713: /* call appropriate ivec function */
714: if (type == GL_MAX)
715: {PCTFS_rvec_max(arg1,arg2,j);}
716: else if (type == GL_MIN)
717: {PCTFS_rvec_min(arg1,arg2,j);}
718: else if (type == GL_MULT)
719: {PCTFS_rvec_mult(arg1,arg2,j);}
720: else if (type == GL_ADD)
721: {PCTFS_rvec_add(arg1,arg2,j);}
722: else if (type == GL_MAX_ABS)
723: {PCTFS_rvec_max_abs(arg1,arg2,j);}
724: else if (type == GL_MIN_ABS)
725: {PCTFS_rvec_min_abs(arg1,arg2,j);}
726: else if (type == GL_EXISTS)
727: {PCTFS_rvec_exists(arg1,arg2,j);}
728: else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_rvec_non_uniform()!");
730: arg1+=j; arg2+=j; i+=j;
731: }
732: return(0);
733: }
735: /***********************************ivec.c*************************************/
736: vfp PCTFS_rvec_fct_addr( PetscInt type)
737: {
738: if (type == NON_UNIFORM)
739: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_non_uniform);}
740: else if (type == GL_MAX)
741: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_max);}
742: else if (type == GL_MIN)
743: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_min);}
744: else if (type == GL_MULT)
745: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_mult);}
746: else if (type == GL_ADD)
747: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_add);}
748: else if (type == GL_MAX_ABS)
749: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_max_abs);}
750: else if (type == GL_MIN_ABS)
751: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_min_abs);}
752: else if (type == GL_EXISTS)
753: {return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_exists);}
755: /* catch all ... not good if we get here */
756: return(NULL);
757: }