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