Actual source code: veccomp.c
slepc-3.7.1 2016-05-27
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: #include <slepc/private/vecimplslepc.h> /*I "slepcvec.h" I*/
24: /* Private MPI datatypes and operators */
25: static MPI_Datatype MPIU_NORM2=0, MPIU_NORM1_AND_2=0;
26: static MPI_Op MPIU_NORM2_SUM=0;
27: static PetscBool VecCompInitialized = PETSC_FALSE;
29: /* Private inline functions */
30: PETSC_STATIC_INLINE void SumNorm2(PetscReal *,PetscReal *,PetscReal *,PetscReal *);
31: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal,PetscReal);
32: PETSC_STATIC_INLINE void AddNorm2(PetscReal *,PetscReal *,PetscReal);
34: #include "veccomp0.h"
37: #include "veccomp0.h"
39: PETSC_STATIC_INLINE void SumNorm2(PetscReal *ssq0,PetscReal *scale0,PetscReal *ssq1,PetscReal *scale1)
40: {
41: PetscReal q;
42: if (*scale0 > *scale1) {
43: q = *scale1/(*scale0);
44: *ssq1 = *ssq0 + q*q*(*ssq1);
45: *scale1 = *scale0;
46: } else {
47: q = *scale0/(*scale1);
48: *ssq1 += q*q*(*ssq0);
49: }
50: }
52: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal ssq,PetscReal scale)
53: {
54: return scale*PetscSqrtReal(ssq);
55: }
57: PETSC_STATIC_INLINE void AddNorm2(PetscReal *ssq,PetscReal *scale,PetscReal x)
58: {
59: PetscReal absx,q;
60: if (x != 0.0) {
61: absx = PetscAbs(x);
62: if (*scale < absx) {
63: q = *scale/absx;
64: *ssq = 1.0 + *ssq*q*q;
65: *scale = absx;
66: } else {
67: q = absx/(*scale);
68: *ssq += q*q;
69: }
70: }
71: }
75: static void SlepcSumNorm2_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
76: {
77: PetscInt i,count = *cnt;
80: if (*datatype == MPIU_NORM2) {
81: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
82: for (i=0; i<count; i++) {
83: SumNorm2(&xin[i*2],&xin[i*2+1],&xout[i*2],&xout[i*2+1]);
84: }
85: } else if (*datatype == MPIU_NORM1_AND_2) {
86: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
87: for (i=0; i<count; i++) {
88: xout[i*3]+= xin[i*3];
89: SumNorm2(&xin[i*3+1],&xin[i*3+2],&xout[i*3+1],&xout[i*3+2]);
90: }
91: } else {
92: (*PetscErrorPrintf)("Can only handle MPIU_NORM* data types");
93: MPI_Abort(MPI_COMM_WORLD,1);
94: }
95: PetscFunctionReturnVoid();
96: }
100: static PetscErrorCode VecNormCompEnd(void)
101: {
105: MPI_Type_free(&MPIU_NORM2);
106: MPI_Type_free(&MPIU_NORM1_AND_2);
107: MPI_Op_free(&MPIU_NORM2_SUM);
108: VecCompInitialized = PETSC_FALSE;
109: return(0);
110: }
114: static PetscErrorCode VecNormCompInit()
115: {
119: MPI_Type_contiguous(sizeof(PetscReal)*2,MPI_BYTE,&MPIU_NORM2);
120: MPI_Type_commit(&MPIU_NORM2);
121: MPI_Type_contiguous(sizeof(PetscReal)*3,MPI_BYTE,&MPIU_NORM1_AND_2);
122: MPI_Type_commit(&MPIU_NORM1_AND_2);
123: MPI_Op_create(SlepcSumNorm2_Local,1,&MPIU_NORM2_SUM);
124: PetscRegisterFinalize(VecNormCompEnd);
125: return(0);
126: }
130: PetscErrorCode VecDestroy_Comp(Vec v)
131: {
132: Vec_Comp *vs = (Vec_Comp*)v->data;
133: PetscInt i;
137: #if defined(PETSC_USE_LOG)
138: PetscLogObjectState((PetscObject)v,"Length=%D",v->map->n);
139: #endif
140: for (i=0;i<vs->nx;i++) {
141: VecDestroy(&vs->x[i]);
142: }
143: if (--vs->n->friends <= 0) {
144: PetscFree(vs->n);
145: }
146: PetscFree(vs->x);
147: PetscFree(vs);
148: return(0);
149: }
151: static struct _VecOps DvOps = {VecDuplicate_Comp, /* 1 */
152: VecDuplicateVecs_Comp,
153: VecDestroyVecs_Comp,
154: VecDot_Comp_MPI,
155: VecMDot_Comp_MPI,
156: VecNorm_Comp_MPI,
157: VecTDot_Comp_MPI,
158: VecMTDot_Comp_MPI,
159: VecScale_Comp,
160: VecCopy_Comp, /* 10 */
161: VecSet_Comp,
162: VecSwap_Comp,
163: VecAXPY_Comp,
164: VecAXPBY_Comp,
165: VecMAXPY_Comp,
166: VecAYPX_Comp,
167: VecWAXPY_Comp,
168: VecAXPBYPCZ_Comp,
169: VecPointwiseMult_Comp,
170: VecPointwiseDivide_Comp,
171: 0, /* 20 */
172: 0,0,
173: 0 /*VecGetArray_Seq*/,
174: VecGetSize_Comp,
175: VecGetLocalSize_Comp,
176: 0/*VecRestoreArray_Seq*/,
177: VecMax_Comp,
178: VecMin_Comp,
179: VecSetRandom_Comp,
180: 0, /* 30 */
181: 0,
182: VecDestroy_Comp,
183: VecView_Comp,
184: 0/*VecPlaceArray_Seq*/,
185: 0/*VecReplaceArray_Seq*/,
186: VecDot_Comp_Seq,
187: 0,
188: VecNorm_Comp_Seq,
189: VecMDot_Comp_Seq,
190: 0, /* 40 */
191: 0,
192: VecReciprocal_Comp,
193: VecConjugate_Comp,
194: 0,0,
195: 0/*VecResetArray_Seq*/,
196: 0,
197: VecMaxPointwiseDivide_Comp,
198: VecPointwiseMax_Comp,
199: VecPointwiseMaxAbs_Comp,
200: VecPointwiseMin_Comp,
201: 0,
202: VecSqrtAbs_Comp,
203: VecAbs_Comp,
204: VecExp_Comp,
205: VecLog_Comp,
206: 0/*VecShift_Comp*/,
207: 0,
208: 0,
209: 0,
210: VecDotNorm2_Comp_MPI
211: };
215: PetscErrorCode VecDuplicateVecs_Comp(Vec w,PetscInt m,Vec *V[])
216: {
218: PetscInt i;
223: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
224: PetscMalloc1(m,V);
225: for (i=0;i<m;i++) { VecDuplicate(w,*V+i); }
226: return(0);
227: }
231: PetscErrorCode VecDestroyVecs_Comp(PetscInt m,Vec v[])
232: {
234: PetscInt i;
238: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
239: for (i=0;i<m;i++) if (v[i]) { VecDestroy(&v[i]); }
240: PetscFree(v);
241: return(0);
242: }
246: static PetscErrorCode VecCreate_Comp_Private(Vec v,Vec *x,PetscInt nx,PetscBool x_to_me,Vec_Comp_N *n)
247: {
248: Vec_Comp *s;
250: PetscInt N=0,lN=0,i,k;
253: if (!VecCompInitialized) {
254: VecCompInitialized = PETSC_TRUE;
255: VecRegister(VECCOMP,VecCreate_Comp);
256: VecNormCompInit();
257: }
259: /* Allocate a new Vec_Comp */
260: if (v->data) { PetscFree(v->data); }
261: PetscNewLog(v,&s);
262: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
263: v->data = (void*)s;
264: v->petscnative = PETSC_FALSE;
266: /* Allocate the array of Vec, if it is needed to be done */
267: if (x_to_me != PETSC_TRUE) {
268: PetscMalloc(sizeof(Vec)*nx,&s->x);
269: PetscMemcpy(s->x,x,sizeof(Vec)*nx);
270: } else s->x = x;
272: s->nx = nx;
274: /* Allocate the shared structure, if it is not given */
275: if (!n) {
276: for (i=0;i<nx;i++) {
277: VecGetSize(x[i],&k);
278: N+= k;
279: VecGetLocalSize(x[i],&k);
280: lN+= k;
281: }
282: PetscNewLog(v,&n);
283: s->n = n;
284: n->n = nx;
285: n->N = N;
286: n->lN = lN;
287: n->friends = 1;
288: } else { /* If not, check in the vector in the shared structure */
289: s->n = n;
290: s->n->friends++;
291: }
293: /* Set the virtual sizes as the real sizes of the vector */
294: VecSetSizes(v,s->n->lN,s->n->N);
296: PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
297: return(0);
298: }
302: PETSC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
303: {
307: VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
308: return(0);
309: }
313: /*@C
314: VecCreateComp - Creates a new vector containing several subvectors,
315: each stored separately.
317: Collective on Vec
319: Input Parameters:
320: + comm - communicator for the new Vec
321: . Nx - array of (initial) global sizes of child vectors
322: . n - number of child vectors
323: . t - type of the child vectors
324: - Vparent - (optional) template vector
326: Output Parameter:
327: . V - new vector
329: Notes:
330: This is similar to PETSc's VecNest but customized for SLEPc's needs. In particular,
331: the number of child vectors can be modified dynamically, with VecCompSetSubVecs().
333: Level: developer
335: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
336: @*/
337: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
338: {
340: Vec *x;
341: PetscInt i;
344: VecCreate(comm,V);
345: PetscMalloc1(n,&x);
346: PetscLogObjectMemory((PetscObject)*V,n*sizeof(Vec));
347: for (i=0;i<n;i++) {
348: VecCreate(comm,&x[i]);
349: VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
350: VecSetType(x[i],t);
351: }
352: VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
353: return(0);
354: }
358: /*@C
359: VecCreateCompWithVecs - Creates a new vector containing several subvectors,
360: each stored separately, from an array of Vecs.
362: Collective on Vec
364: Input Parameters:
365: + x - array of Vecs
366: . n - number of child vectors
367: - Vparent - (optional) template vector
369: Output Parameter:
370: . V - new vector
372: Level: developer
374: .seealso: VecCreateComp()
375: @*/
376: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
377: {
379: PetscInt i;
382: VecCreate(PetscObjectComm((PetscObject)x[0]),V);
383: for (i=0;i<n;i++) {
384: PetscObjectReference((PetscObject)x[i]);
385: }
386: VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
387: return(0);
388: }
392: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
393: {
395: Vec *x;
396: PetscInt i;
397: Vec_Comp *s = (Vec_Comp*)win->data;
400: SlepcValidVecComp(win);
401: VecCreate(PetscObjectComm((PetscObject)win),V);
402: PetscMalloc1(s->nx,&x);
403: PetscLogObjectMemory((PetscObject)*V,s->nx*sizeof(Vec));
404: for (i=0;i<s->nx;i++) {
405: if (s->x[i]) {
406: VecDuplicate(s->x[i],&x[i]);
407: } else x[i] = NULL;
408: }
409: VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
410: return(0);
411: }
415: /*@C
416: VecCompGetSubVecs - Returns the entire array of vectors defining a
417: compound vector.
419: Collective on Vec
421: Input Parameter:
422: . win - compound vector
424: Output Parameters:
425: + n - number of child vectors
426: - x - array of child vectors
428: Level: developer
430: .seealso: VecCreateComp()
431: @*/
432: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
433: {
434: Vec_Comp *s = (Vec_Comp*)win->data;
437: SlepcValidVecComp(win);
438: if (x) *x = s->x;
439: if (n) *n = s->nx;
440: return(0);
441: }
445: /*@C
446: VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
447: of replaces the subvectors.
449: Collective on Vec
451: Input Parameters:
452: + win - compound vector
453: . n - number of child vectors
454: - x - array of child vectors
456: Level: developer
458: .seealso: VecCreateComp()
459: @*/
460: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
461: {
462: Vec_Comp *s = (Vec_Comp*)win->data;
466: if (x) {
467: if (n > s->nx) {
468: PetscFree(s->x);
469: PetscMalloc(sizeof(Vec)*n,&s->x);
470: }
471: PetscMemcpy(s->x,x,sizeof(Vec)*n);
472: s->nx = n;
473: }
474: s->n->n = n;
475: return(0);
476: }
480: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
481: {
483: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
484: PetscInt i;
487: SlepcValidVecComp(v);
488: SlepcValidVecComp(w);
489: for (i=0;i<vs->n->n;i++) {
490: VecAXPY(vs->x[i],alpha,ws->x[i]);
491: }
492: return(0);
493: }
497: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
498: {
500: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
501: PetscInt i;
504: SlepcValidVecComp(v);
505: SlepcValidVecComp(w);
506: for (i=0;i<vs->n->n;i++) {
507: VecAYPX(vs->x[i],alpha,ws->x[i]);
508: }
509: return(0);
510: }
514: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
515: {
517: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
518: PetscInt i;
521: SlepcValidVecComp(v);
522: SlepcValidVecComp(w);
523: for (i=0;i<vs->n->n;i++) {
524: VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
525: }
526: return(0);
527: }
531: PetscErrorCode VecMAXPY_Comp(Vec v,PetscInt n,const PetscScalar *alpha,Vec *w)
532: {
534: Vec_Comp *vs = (Vec_Comp*)v->data;
535: Vec *wx;
536: PetscInt i,j;
539: SlepcValidVecComp(v);
540: for (i=0;i<n;i++) SlepcValidVecComp(w[i]);
542: PetscMalloc(sizeof(Vec)*n,&wx);
544: for (j=0;j<vs->n->n;j++) {
545: for (i=0;i<n;i++) wx[i] = ((Vec_Comp*)w[i]->data)->x[j];
546: VecMAXPY(vs->x[j],n,alpha,wx);
547: }
549: PetscFree(wx);
550: return(0);
551: }
555: PetscErrorCode VecWAXPY_Comp(Vec v,PetscScalar alpha,Vec w,Vec z)
556: {
558: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
559: PetscInt i;
562: SlepcValidVecComp(v);
563: SlepcValidVecComp(w);
564: SlepcValidVecComp(z);
565: for (i=0;i<vs->n->n;i++) {
566: VecWAXPY(vs->x[i],alpha,ws->x[i],zs->x[i]);
567: }
568: return(0);
569: }
573: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
574: {
575: PetscErrorCode ierr;
576: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
577: PetscInt i;
580: SlepcValidVecComp(v);
581: SlepcValidVecComp(w);
582: SlepcValidVecComp(z);
583: for (i=0;i<vs->n->n;i++) {
584: VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
585: }
586: return(0);
587: }
591: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
592: {
593: Vec_Comp *vs = (Vec_Comp*)v->data;
596: SlepcValidVecComp(v);
598: *size = vs->n->N;
599: return(0);
600: }
604: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
605: {
606: Vec_Comp *vs = (Vec_Comp*)v->data;
609: SlepcValidVecComp(v);
611: *size = vs->n->lN;
612: return(0);
613: }
617: PetscErrorCode VecMax_Comp(Vec v,PetscInt *idx,PetscReal *z)
618: {
620: Vec_Comp *vs = (Vec_Comp*)v->data;
621: PetscInt idxp,s=0,s0;
622: PetscReal zp,z0;
623: PetscInt i;
626: SlepcValidVecComp(v);
627: if (!idx && !z) return(0);
629: if (vs->n->n > 0) {
630: VecMax(vs->x[0],idx?&idxp:NULL,&zp);
631: } else {
632: zp = PETSC_MIN_REAL;
633: if (idx) idxp = -1;
634: }
635: for (i=1;i<vs->n->n;i++) {
636: VecGetSize(vs->x[i-1],&s0);
637: s += s0;
638: VecMax(vs->x[i],idx?&idxp:NULL,&z0);
639: if (zp < z0) {
640: if (idx) *idx = s+idxp;
641: zp = z0;
642: }
643: }
644: if (z) *z = zp;
645: return(0);
646: }
650: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
651: {
653: Vec_Comp *vs = (Vec_Comp*)v->data;
654: PetscInt idxp,s=0,s0;
655: PetscReal zp,z0;
656: PetscInt i;
659: SlepcValidVecComp(v);
660: if (!idx && !z) return(0);
662: if (vs->n->n > 0) {
663: VecMin(vs->x[0],idx?&idxp:NULL,&zp);
664: } else {
665: zp = PETSC_MAX_REAL;
666: if (idx) idxp = -1;
667: }
668: for (i=1;i<vs->n->n;i++) {
669: VecGetSize(vs->x[i-1],&s0);
670: s += s0;
671: VecMin(vs->x[i],idx?&idxp:NULL,&z0);
672: if (zp > z0) {
673: if (idx) *idx = s+idxp;
674: zp = z0;
675: }
676: }
677: if (z) *z = zp;
678: return(0);
679: }
683: PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v,Vec w,PetscReal *m)
684: {
686: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
687: PetscReal work;
688: PetscInt i;
691: SlepcValidVecComp(v);
692: SlepcValidVecComp(w);
693: if (!m || vs->n->n == 0) return(0);
694: VecMaxPointwiseDivide(vs->x[0],ws->x[0],m);
695: for (i=1;i<vs->n->n;i++) {
696: VecMaxPointwiseDivide(vs->x[i],ws->x[i],&work);
697: *m = PetscMax(*m,work);
698: }
699: return(0);
700: }
708: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v) \
709: { \
710: PetscErrorCode ierr; \
711: Vec_Comp *vs = (Vec_Comp*)v->data; \
712: PetscInt i; \
713: \
715: SlepcValidVecComp(v); \
716: for (i=0;i<vs->n->n;i++) { \
717: __COMPOSE2__(Vec,NAME)(vs->x[i]); \
718: } \
719: return(0);\
720: }
724: __FUNC_TEMPLATE1__(Conjugate)
728: __FUNC_TEMPLATE1__(Reciprocal)
732: __FUNC_TEMPLATE1__(SqrtAbs)
736: __FUNC_TEMPLATE1__(Abs)
740: __FUNC_TEMPLATE1__(Exp)
744: __FUNC_TEMPLATE1__(Log)
748: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,T0 __a) \
749: { \
750: PetscErrorCode ierr; \
751: Vec_Comp *vs = (Vec_Comp*)v->data; \
752: PetscInt i; \
753: \
755: SlepcValidVecComp(v); \
756: for (i=0;i<vs->n->n;i++) { \
757: __COMPOSE2__(Vec,NAME)(vs->x[i],__a); \
758: } \
759: return(0);\
760: }
764: __FUNC_TEMPLATE2__(Set,PetscScalar)
768: __FUNC_TEMPLATE2__(View,PetscViewer)
772: __FUNC_TEMPLATE2__(Scale,PetscScalar)
776: __FUNC_TEMPLATE2__(SetRandom,PetscRandom)
780: __FUNC_TEMPLATE2__(Shift,PetscScalar)
784: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w) \
785: { \
786: PetscErrorCode ierr; \
787: Vec_Comp *vs = (Vec_Comp*)v->data,\
788: *ws = (Vec_Comp*)w->data; \
789: PetscInt i; \
790: \
792: SlepcValidVecComp(v); \
793: SlepcValidVecComp(w); \
794: for (i=0;i<vs->n->n;i++) { \
795: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i]); \
796: } \
797: return(0);\
798: }
802: __FUNC_TEMPLATE3__(Copy)
806: __FUNC_TEMPLATE3__(Swap)
810: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w,Vec z) \
811: { \
812: PetscErrorCode ierr; \
813: Vec_Comp *vs = (Vec_Comp*)v->data, \
814: *ws = (Vec_Comp*)w->data, \
815: *zs = (Vec_Comp*)z->data; \
816: PetscInt i; \
817: \
819: SlepcValidVecComp(v); \
820: SlepcValidVecComp(w); \
821: SlepcValidVecComp(z); \
822: for (i=0;i<vs->n->n;i++) { \
823: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i],zs->x[i]); \
824: } \
825: return(0);\
826: }
830: __FUNC_TEMPLATE4__(PointwiseMax)
834: __FUNC_TEMPLATE4__(PointwiseMaxAbs)
838: __FUNC_TEMPLATE4__(PointwiseMin)
842: __FUNC_TEMPLATE4__(PointwiseMult)
846: __FUNC_TEMPLATE4__(PointwiseDivide)