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