xref: /petsc/src/ksp/pc/impls/tfs/ivec.c (revision 009bbdc485cd9ad46be9940d3549e2dde9cdc322)
1 
2 /**********************************ivec.c**************************************
3 
4 Author: Henry M. Tufo III
5 
6 e-mail: hmt@cs.brown.edu
7 
8 snail-mail:
9 Division of Applied Mathematics
10 Brown University
11 Providence, RI 02912
12 
13 Last Modification:
14 6.21.97
15 ***********************************ivec.c*************************************/
16 
17 #include <../src/ksp/pc/impls/tfs/tfs.h>
18 
19 /* sorting args ivec.c ivec.c ... */
20 #define   SORT_OPT     6
21 #define   SORT_STACK   50000
22 
23 
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];
27 
28 /***********************************ivec.c*************************************/
29 PetscInt *PCTFS_ivec_copy(PetscInt *arg1, PetscInt *arg2, PetscInt n)
30 {
31   while (n--) { *arg1++ = *arg2++; }
32   return(arg1);
33 }
34 
35 /***********************************ivec.c*************************************/
36 PetscErrorCode PCTFS_ivec_zero(PetscInt *arg1, PetscInt n)
37 {
38   PetscFunctionBegin;
39   while (n--) { *arg1++ = 0; }
40   PetscFunctionReturn(0);
41 }
42 
43 /***********************************ivec.c*************************************/
44 PetscErrorCode PCTFS_ivec_set(PetscInt *arg1, PetscInt arg2, PetscInt n)
45 {
46   PetscFunctionBegin;
47   while (n--) { *arg1++ = arg2; }
48   PetscFunctionReturn(0);
49 }
50 
51 /***********************************ivec.c*************************************/
52 PetscErrorCode PCTFS_ivec_max(PetscInt *arg1, PetscInt *arg2, PetscInt n)
53 {
54   PetscFunctionBegin;
55   while (n--) { *arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++; }
56   PetscFunctionReturn(0);
57 }
58 
59 /***********************************ivec.c*************************************/
60 PetscErrorCode PCTFS_ivec_min(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
61 {
62   PetscFunctionBegin;
63   while (n--) { *(arg1) = PetscMin(*arg1,*arg2); arg1++; arg2++; }
64   PetscFunctionReturn(0);
65 }
66 
67 /***********************************ivec.c*************************************/
68 PetscErrorCode PCTFS_ivec_mult(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
69 {
70   PetscFunctionBegin;
71   while (n--) { *arg1++ *= *arg2++; }
72   PetscFunctionReturn(0);
73 }
74 
75 /***********************************ivec.c*************************************/
76 PetscErrorCode PCTFS_ivec_add(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
77 {
78   PetscFunctionBegin;
79   while (n--) { *arg1++ += *arg2++; }
80   PetscFunctionReturn(0);
81 }
82 
83 /***********************************ivec.c*************************************/
84 PetscErrorCode PCTFS_ivec_lxor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
85 {
86   PetscFunctionBegin;
87   while (n--) { *arg1=((*arg1 || *arg2) && !(*arg1 && *arg2)); arg1++; arg2++; }
88   PetscFunctionReturn(0);
89 }
90 
91 /***********************************ivec.c*************************************/
92 PetscErrorCode PCTFS_ivec_xor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
93 {
94   PetscFunctionBegin;
95   while (n--) { *arg1++ ^= *arg2++; }
96   PetscFunctionReturn(0);
97 }
98 
99 /***********************************ivec.c*************************************/
100 PetscErrorCode PCTFS_ivec_or(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
101 {
102   PetscFunctionBegin;
103   while (n--) { *arg1++ |= *arg2++; }
104   PetscFunctionReturn(0);
105 }
106 
107 /***********************************ivec.c*************************************/
108 PetscErrorCode PCTFS_ivec_lor(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
109 {
110   PetscFunctionBegin;
111   while (n--) { *arg1 = (*arg1 || *arg2); arg1++; arg2++; }
112   PetscFunctionReturn(0);
113 }
114 
115 /***********************************ivec.c*************************************/
116 PetscErrorCode PCTFS_ivec_and(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
117 {
118   PetscFunctionBegin;
119   while (n--) { *arg1++ &= *arg2++; }
120   PetscFunctionReturn(0);
121 }
122 
123 /***********************************ivec.c*************************************/
124 PetscErrorCode PCTFS_ivec_land(PetscInt *arg1,  PetscInt *arg2,  PetscInt n)
125 {
126   PetscFunctionBegin;
127   while (n--) { *arg1 = (*arg1 && *arg2); arg1++; arg2++; }
128   PetscFunctionReturn(0);
129 }
130 
131 /***********************************ivec.c*************************************/
132 PetscErrorCode PCTFS_ivec_and3(PetscInt *arg1,  PetscInt *arg2,  PetscInt *arg3, PetscInt n)
133 {
134   PetscFunctionBegin;
135   while (n--) { *arg1++ = (*arg2++ & *arg3++); }
136   PetscFunctionReturn(0);
137 }
138 
139 /***********************************ivec.c*************************************/
140 PetscInt PCTFS_ivec_sum(PetscInt *arg1,  PetscInt n)
141 {
142   PetscInt tmp = 0;
143 
144   while (n--) {tmp += *arg1++;}
145   return(tmp);
146 }
147 
148 /***********************************ivec.c*************************************/
149 PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *arg1, PetscInt *arg2,  PetscInt n,  PetscInt *arg3)
150 {
151   PetscInt i, j, type;
152 
153   PetscFunctionBegin;
154   /* LATER: if we're really motivated we can sort and then unsort */
155   for (i=0;i<n;) {
156     /* clump 'em for now */
157     j=i+1;
158     type = arg3[i];
159     while ((j<n)&&(arg3[j]==type)) { j++; }
160 
161     /* how many together */
162     j -= i;
163 
164     /* call appropriate ivec function */
165     if (type == GL_MAX)        { PCTFS_ivec_max(arg1,arg2,j); }
166     else if (type == GL_MIN)   { PCTFS_ivec_min(arg1,arg2,j); }
167     else if (type == GL_MULT)  { PCTFS_ivec_mult(arg1,arg2,j); }
168     else if (type == GL_ADD)   { PCTFS_ivec_add(arg1,arg2,j); }
169     else if (type == GL_B_XOR) { PCTFS_ivec_xor(arg1,arg2,j); }
170     else if (type == GL_B_OR)  { PCTFS_ivec_or(arg1,arg2,j); }
171     else if (type == GL_B_AND) { PCTFS_ivec_and(arg1,arg2,j); }
172     else if (type == GL_L_XOR) { PCTFS_ivec_lxor(arg1,arg2,j); }
173     else if (type == GL_L_OR)  { PCTFS_ivec_lor(arg1,arg2,j); }
174     else if (type == GL_L_AND) { PCTFS_ivec_land(arg1,arg2,j); }
175     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_ivec_non_uniform()!");
176 
177     arg1+=j; arg2+=j; i+=j;
178   }
179   PetscFunctionReturn(0);
180 }
181 
182 /***********************************ivec.c*************************************/
183 vfp PCTFS_ivec_fct_addr(PetscInt type)
184 {
185   if (type == NON_UNIFORM)   { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_non_uniform); }
186   else if (type == GL_MAX)   { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_max); }
187   else if (type == GL_MIN)   { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_min); }
188   else if (type == GL_MULT)  { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_mult); }
189   else if (type == GL_ADD)   { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_add); }
190   else if (type == GL_B_XOR) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_xor); }
191   else if (type == GL_B_OR)  { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_or); }
192   else if (type == GL_B_AND) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_and); }
193   else if (type == GL_L_XOR) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_lxor); }
194   else if (type == GL_L_OR)  { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_lor); }
195   else if (type == GL_L_AND) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_ivec_land); }
196 
197   /* catch all ... not good if we get here */
198   return(NULL);
199 }
200 
201 /******************************************************************************/
202 PetscErrorCode PCTFS_ivec_sort(PetscInt *ar,  PetscInt size)
203 {
204   PetscInt *pi, *pj, temp;
205   PetscInt **top_a = (PetscInt **) offset_stack;
206   PetscInt *top_s = size_stack, *bottom_s = size_stack;
207 
208 
209   /* we're really interested in the offset of the last element */
210   /* ==> length of the list is now size + 1                    */
211   size--;
212 
213   /* do until we're done ... return when stack is exhausted */
214   for (;;) {
215     /* if list is large enough use quicksort partition exchange code */
216     if (size > SORT_OPT) {
217       /* start up pointer at element 1 and down at size     */
218       pi = ar+1;
219       pj = ar+size;
220 
221       /* find middle element in list and swap w/ element 1 */
222       SWAP(*(ar+(size>>1)),*pi)
223 
224       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
225       /* note ==> pivot_value in index 0                   */
226       if (*pi > *pj) { SWAP(*pi,*pj) }
227       if (*ar > *pj) { SWAP(*ar,*pj) }
228       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) }
229 
230       /* partition about pivot_value ...                              */
231       /* note lists of length 2 are not guaranteed to be sorted */
232       for (;;) {
233         /* walk up ... and down ... swap if equal to pivot! */
234         do pi++; while (*pi<*ar);
235         do pj--; while (*pj>*ar);
236 
237         /* if we've crossed we're done */
238         if (pj<pi) break;
239 
240         /* else swap */
241         SWAP(*pi,*pj)
242       }
243 
244       /* place pivot_value in it's correct location */
245       SWAP(*ar,*pj)
246 
247       /* test stack_size to see if we've exhausted our stack */
248       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort() :: STACK EXHAUSTED!!!");
249 
250       /* push right hand child iff length > 1 */
251       if ((*top_s = size-((PetscInt) (pi-ar)))) {
252         *(top_a++) = pi;
253         size -= *top_s+2;
254         top_s++;
255       } else if (size -= *top_s+2) {;} /* set up for next loop iff there is something to do */
256       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
257         ar = *(--top_a);
258         size = *(--top_s);
259       }
260     } else { /* else sort small list directly then pop another off stack */
261 
262       /* insertion sort for bottom */
263       for (pj=ar+1;pj<=ar+size;pj++) {
264         temp = *pj;
265         for (pi=pj-1;pi>=ar;pi--) {
266           if (*pi <= temp) break;
267           *(pi+1)=*pi;
268         }
269         *(pi+1)=temp;
270       }
271 
272       /* check to see if stack is exhausted ==> DONE */
273       if (top_s==bottom_s) PetscFunctionReturn(0);
274 
275       /* else pop another list from the stack */
276       ar = *(--top_a);
277       size = *(--top_s);
278     }
279   }
280   PetscFunctionReturn(0);
281 }
282 
283 /******************************************************************************/
284 PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *ar,  PetscInt *ar2,  PetscInt size)
285 {
286   PetscInt *pi, *pj, temp, temp2;
287   PetscInt **top_a = (PetscInt **)offset_stack;
288   PetscInt *top_s = size_stack, *bottom_s = size_stack;
289   PetscInt *pi2, *pj2;
290   PetscInt mid;
291 
292   PetscFunctionBegin;
293   /* we're really interested in the offset of the last element */
294   /* ==> length of the list is now size + 1                    */
295   size--;
296 
297   /* do until we're done ... return when stack is exhausted */
298   for (;;) {
299 
300     /* if list is large enough use quicksort partition exchange code */
301     if (size > SORT_OPT) {
302 
303       /* start up pointer at element 1 and down at size     */
304       mid = size>>1;
305       pi = ar+1;
306       pj = ar+mid;
307       pi2 = ar2+1;
308       pj2 = ar2+mid;
309 
310       /* find middle element in list and swap w/ element 1 */
311       SWAP(*pi,*pj)
312       SWAP(*pi2,*pj2)
313 
314       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
315       /* note ==> pivot_value in index 0                   */
316       pj = ar+size;
317       pj2 = ar2+size;
318       if (*pi > *pj) { SWAP(*pi,*pj) SWAP(*pi2,*pj2) }
319       if (*ar > *pj) { SWAP(*ar,*pj) SWAP(*ar2,*pj2) }
320       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) SWAP(*(ar2),*(ar2+1)) }
321 
322       /* partition about pivot_value ...                              */
323       /* note lists of length 2 are not guaranteed to be sorted */
324       for (;;) {
325         /* walk up ... and down ... swap if equal to pivot! */
326         do { pi++; pi2++; } while (*pi<*ar);
327         do { pj--; pj2--; } while (*pj>*ar);
328 
329         /* if we've crossed we're done */
330         if (pj<pi) break;
331 
332         /* else swap */
333         SWAP(*pi,*pj)
334         SWAP(*pi2,*pj2)
335       }
336 
337       /* place pivot_value in it's correct location */
338       SWAP(*ar,*pj)
339       SWAP(*ar2,*pj2)
340 
341       /* test stack_size to see if we've exhausted our stack */
342       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion() :: STACK EXHAUSTED!!!");
343 
344       /* push right hand child iff length > 1 */
345       if ((*top_s = size-((PetscInt) (pi-ar)))) {
346         *(top_a++) = pi;
347         *(top_a++) = pi2;
348         size -= *top_s+2;
349         top_s++;
350       } else if (size -= *top_s+2) {;} /* set up for next loop iff there is something to do */
351       else  { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
352         ar2 = *(--top_a);
353         ar  = *(--top_a);
354         size = *(--top_s);
355       }
356     } else { /* else sort small list directly then pop another off stack */
357 
358       /* insertion sort for bottom */
359       for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++) {
360         temp = *pj;
361         temp2 = *pj2;
362         for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--) {
363           if (*pi <= temp) break;
364           *(pi+1)=*pi;
365           *(pi2+1)=*pi2;
366         }
367         *(pi+1)=temp;
368         *(pi2+1)=temp2;
369       }
370 
371       /* check to see if stack is exhausted ==> DONE */
372       if (top_s==bottom_s) PetscFunctionReturn(0);
373 
374       /* else pop another list from the stack */
375       ar2 = *(--top_a);
376       ar  = *(--top_a);
377       size = *(--top_s);
378     }
379   }
380   PetscFunctionReturn(0);
381 }
382 
383 /******************************************************************************/
384 PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *ar,  PetscInt **ar2, PetscInt size)
385 {
386   PetscInt *pi, *pj, temp, *ptr;
387   PetscInt **top_a = (PetscInt **)offset_stack;
388   PetscInt *top_s = size_stack, *bottom_s = size_stack;
389   PetscInt **pi2, **pj2;
390   PetscInt mid;
391 
392   PetscFunctionBegin;
393   /* we're really interested in the offset of the last element */
394   /* ==> length of the list is now size + 1                    */
395   size--;
396 
397   /* do until we're done ... return when stack is exhausted */
398   for (;;) {
399 
400     /* if list is large enough use quicksort partition exchange code */
401     if (size > SORT_OPT) {
402 
403       /* start up pointer at element 1 and down at size     */
404       mid = size>>1;
405       pi = ar+1;
406       pj = ar+mid;
407       pi2 = ar2+1;
408       pj2 = ar2+mid;
409 
410       /* find middle element in list and swap w/ element 1 */
411       SWAP(*pi,*pj)
412       P_SWAP(*pi2,*pj2)
413 
414       /* order element 0,1,size-1 st {M,L,...,U} w/L<=M<=U */
415       /* note ==> pivot_value in index 0                   */
416       pj = ar+size;
417       pj2 = ar2+size;
418       if (*pi > *pj) { SWAP(*pi,*pj) P_SWAP(*pi2,*pj2) }
419       if (*ar > *pj) { SWAP(*ar,*pj) P_SWAP(*ar2,*pj2) }
420       else if (*pi > *ar) { SWAP(*(ar),*(ar+1)) P_SWAP(*(ar2),*(ar2+1)) }
421 
422       /* partition about pivot_value ...                              */
423       /* note lists of length 2 are not guaranteed to be sorted */
424       for (;;) {
425 
426         /* walk up ... and down ... swap if equal to pivot! */
427         do {pi++; pi2++;} while (*pi<*ar);
428         do {pj--; pj2--;} while (*pj>*ar);
429 
430         /* if we've crossed we're done */
431         if (pj<pi) break;
432 
433         /* else swap */
434         SWAP(*pi,*pj)
435         P_SWAP(*pi2,*pj2)
436       }
437 
438       /* place pivot_value in it's correct location */
439       SWAP(*ar,*pj)
440       P_SWAP(*ar2,*pj2)
441 
442       /* test stack_size to see if we've exhausted our stack */
443       if (top_s-bottom_s >= SORT_STACK) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_ivec_sort_companion_hack() :: STACK EXHAUSTED!!!");
444 
445       /* push right hand child iff length > 1 */
446       if ((*top_s = size-((PetscInt) (pi-ar)))) {
447         *(top_a++) = pi;
448         *(top_a++) = (PetscInt*) pi2;
449         size -= *top_s+2;
450         top_s++;
451       } else if (size -= *top_s+2) {;} /* set up for next loop iff there is something to do */
452       else { /* might as well pop - note NR_OPT >=2 ==> we're ok! */
453         ar2 = (PetscInt **) *(--top_a);
454         ar  = *(--top_a);
455         size = *(--top_s);
456       }
457     }
458 
459 
460     else { /* else sort small list directly then pop another off stack */
461       /* insertion sort for bottom */
462       for (pj=ar+1, pj2=ar2+1;pj<=ar+size;pj++,pj2++) {
463         temp = *pj;
464         ptr = *pj2;
465         for (pi=pj-1,pi2=pj2-1;pi>=ar;pi--,pi2--) {
466           if (*pi <= temp) break;
467           *(pi+1)=*pi;
468           *(pi2+1)=*pi2;
469         }
470         *(pi+1)=temp;
471         *(pi2+1)=ptr;
472       }
473 
474       /* check to see if stack is exhausted ==> DONE */
475       if (top_s==bottom_s) PetscFunctionReturn(0);
476 
477       /* else pop another list from the stack */
478       ar2 = (PetscInt **)*(--top_a);
479       ar  = *(--top_a);
480       size = *(--top_s);
481     }
482   }
483   PetscFunctionReturn(0);
484 }
485 
486 /******************************************************************************/
487 PetscErrorCode PCTFS_SMI_sort(void *ar1, void *ar2, PetscInt size, PetscInt type)
488 {
489   PetscFunctionBegin;
490   if (type == SORT_INTEGER) {
491     if (ar2) PCTFS_ivec_sort_companion((PetscInt*)ar1,(PetscInt*)ar2,size);
492     else PCTFS_ivec_sort((PetscInt*)ar1,size);
493   } else if (type == SORT_INT_PTR) {
494     if (ar2) PCTFS_ivec_sort_companion_hack((PetscInt*)ar1,(PetscInt **)ar2,size);
495     else PCTFS_ivec_sort((PetscInt*)ar1,size);
496   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCTFS_SMI_sort only does SORT_INTEGER!");
497   PetscFunctionReturn(0);
498 }
499 
500 /***********************************ivec.c*************************************/
501 PetscInt PCTFS_ivec_linear_search(PetscInt item,  PetscInt *list,  PetscInt n)
502 {
503   PetscInt tmp = n-1;
504 
505   PetscFunctionBegin;
506   while (n--)  { if (*list++ == item) { return(tmp-n); } }
507   return(-1);
508 }
509 
510 /***********************************ivec.c*************************************/
511 PetscInt PCTFS_ivec_binary_search(PetscInt item,  PetscInt *list,  PetscInt rh)
512 {
513   PetscInt mid, lh=0;
514 
515   rh--;
516   while (lh<=rh) {
517     mid = (lh+rh)>>1;
518     if (*(list+mid) == item) { return(mid); }
519     if (*(list+mid) > item) { rh = mid-1; }
520     else { lh = mid+1; }
521   }
522   return(-1);
523 }
524 
525 /*********************************ivec.c*************************************/
526 PetscErrorCode PCTFS_rvec_copy(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
527 {
528   PetscFunctionBegin;
529   while (n--) { *arg1++ = *arg2++; }
530   PetscFunctionReturn(0);
531 }
532 
533 /*********************************ivec.c*************************************/
534 PetscErrorCode PCTFS_rvec_zero(PetscScalar *arg1,  PetscInt n)
535 {
536   PetscFunctionBegin;
537   while (n--) { *arg1++ = 0.0; }
538   PetscFunctionReturn(0);
539 }
540 
541 /***********************************ivec.c*************************************/
542 PetscErrorCode PCTFS_rvec_one(PetscScalar *arg1,  PetscInt n)
543 {
544   PetscFunctionBegin;
545   while (n--) { *arg1++ = 1.0; }
546   PetscFunctionReturn(0);
547 }
548 
549 /***********************************ivec.c*************************************/
550 PetscErrorCode PCTFS_rvec_set(PetscScalar *arg1,  PetscScalar arg2,  PetscInt n)
551 {
552   PetscFunctionBegin;
553   while (n--) { *arg1++ = arg2; }
554   PetscFunctionReturn(0);
555 }
556 
557 /***********************************ivec.c*************************************/
558 PetscErrorCode PCTFS_rvec_scale(PetscScalar *arg1,  PetscScalar arg2,  PetscInt n)
559 {
560   PetscFunctionBegin;
561   while (n--) { *arg1++ *= arg2; }
562   PetscFunctionReturn(0);
563 }
564 
565 /*********************************ivec.c*************************************/
566 PetscErrorCode PCTFS_rvec_add(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
567 {
568   PetscFunctionBegin;
569   while (n--) { *arg1++ += *arg2++; }
570   PetscFunctionReturn(0);
571 }
572 
573 /*********************************ivec.c*************************************/
574 PetscErrorCode PCTFS_rvec_mult(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
575 {
576   PetscFunctionBegin;
577   while (n--) { *arg1++ *= *arg2++; }
578   PetscFunctionReturn(0);
579 }
580 
581 /*********************************ivec.c*************************************/
582 PetscErrorCode PCTFS_rvec_max(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
583 {
584   PetscFunctionBegin;
585   while (n--) { *arg1 = PetscMax(*arg1,*arg2); arg1++; arg2++; }
586   PetscFunctionReturn(0);
587 }
588 
589 /*********************************ivec.c*************************************/
590 PetscErrorCode PCTFS_rvec_max_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
591 {
592   PetscFunctionBegin;
593   while (n--) { *arg1 = MAX_FABS(*arg1,*arg2); arg1++; arg2++; }
594   PetscFunctionReturn(0);
595 }
596 
597 /*********************************ivec.c*************************************/
598 PetscErrorCode PCTFS_rvec_min(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
599 {
600   PetscFunctionBegin;
601   while (n--) { *arg1 = PetscMin(*arg1,*arg2); arg1++; arg2++; }
602   PetscFunctionReturn(0);
603 }
604 
605 /*********************************ivec.c*************************************/
606 PetscErrorCode PCTFS_rvec_min_abs(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
607 {
608   PetscFunctionBegin;
609   while (n--) { *arg1 = MIN_FABS(*arg1,*arg2); arg1++; arg2++; }
610   PetscFunctionReturn(0);
611 }
612 
613 /*********************************ivec.c*************************************/
614 PetscErrorCode PCTFS_rvec_exists(PetscScalar *arg1,  PetscScalar *arg2,  PetscInt n)
615 {
616   PetscFunctionBegin;
617   while (n--) { *arg1 = EXISTS(*arg1,*arg2); arg1++; arg2++; }
618   PetscFunctionReturn(0);
619 }
620 
621 /***********************************ivec.c*************************************/
622 PetscErrorCode PCTFS_rvec_non_uniform(PetscScalar *arg1, PetscScalar *arg2,  PetscInt n,  PetscInt *arg3)
623 {
624   PetscInt i, j, type;
625 
626   PetscFunctionBegin;
627   /* LATER: if we're really motivated we can sort and then unsort */
628   for (i=0;i<n;) {
629 
630     /* clump 'em for now */
631     j=i+1;
632     type = arg3[i];
633     while ((j<n)&&(arg3[j]==type)) { j++; }
634 
635     /* how many together */
636     j -= i;
637 
638     /* call appropriate ivec function */
639     if (type == GL_MAX)          { PCTFS_rvec_max(arg1,arg2,j); }
640     else if (type == GL_MIN)     { PCTFS_rvec_min(arg1,arg2,j); }
641     else if (type == GL_MULT)    { PCTFS_rvec_mult(arg1,arg2,j); }
642     else if (type == GL_ADD)     { PCTFS_rvec_add(arg1,arg2,j); }
643     else if (type == GL_MAX_ABS) { PCTFS_rvec_max_abs(arg1,arg2,j); }
644     else if (type == GL_MIN_ABS) { PCTFS_rvec_min_abs(arg1,arg2,j); }
645     else if (type == GL_EXISTS)  { PCTFS_rvec_exists(arg1,arg2,j); }
646     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"unrecognized type passed to PCTFS_rvec_non_uniform()!");
647 
648     arg1+=j; arg2+=j; i+=j;
649   }
650   PetscFunctionReturn(0);
651 }
652 
653 /***********************************ivec.c*************************************/
654 vfp PCTFS_rvec_fct_addr(PetscInt type)
655 {
656   if (type == NON_UNIFORM)     { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_non_uniform); }
657   else if (type == GL_MAX)     { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_max); }
658   else if (type == GL_MIN)     { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_min); }
659   else if (type == GL_MULT)    { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_mult); }
660   else if (type == GL_ADD)     { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_add); }
661   else if (type == GL_MAX_ABS) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_max_abs); }
662   else if (type == GL_MIN_ABS) { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_min_abs); }
663   else if (type == GL_EXISTS)  { return((PetscErrorCode (*)(void*, void *, PetscInt, ...))&PCTFS_rvec_exists); }
664 
665   /* catch all ... not good if we get here */
666   return(NULL);
667 }
668 
669 
670 
671 
672 
673 
674