Actual source code: veccomp.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2013, 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: #if !defined(PETSC_USE_DYNAMIC_LIBRARIES)
28: static PetscBool VecCompInitialized = PETSC_FALSE;
29: #endif
31: /* Private inline functions */
32: PETSC_STATIC_INLINE void SumNorm2(PetscReal *,PetscReal *,PetscReal *,PetscReal *);
33: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal,PetscReal);
34: PETSC_STATIC_INLINE void AddNorm2(PetscReal *,PetscReal *,PetscReal);
36: #include "veccomp0.h"
39: #include "veccomp0.h"
41: PETSC_STATIC_INLINE void SumNorm2(PetscReal *ssq0,PetscReal *scale0,PetscReal *ssq1,PetscReal *scale1)
42: {
43: PetscReal q;
44: if (*scale0 > *scale1) {
45: q = *scale1/(*scale0);
46: *ssq1 = *ssq0 + q*q*(*ssq1);
47: *scale1 = *scale0;
48: } else {
49: q = *scale0/(*scale1);
50: *ssq1 += q*q*(*ssq0);
51: }
52: }
54: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal ssq,PetscReal scale)
55: {
56: return scale*PetscSqrtReal(ssq);
57: }
59: PETSC_STATIC_INLINE void AddNorm2(PetscReal *ssq,PetscReal *scale,PetscReal x)
60: {
61: PetscReal absx,q;
62: if (x != 0.0) {
63: absx = PetscAbs(x);
64: if (*scale < absx) {
65: q = *scale/absx;
66: *ssq = 1.0 + *ssq*q*q;
67: *scale = absx;
68: } else {
69: q = absx/(*scale);
70: *ssq += q*q;
71: }
72: }
73: }
75: #if !defined(PETSC_USE_DYNAMIC_LIBRARIES)
79: static void SlepcSumNorm2_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
80: {
81: PetscInt i,count = *cnt;
84: if (*datatype == MPIU_NORM2) {
85: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
86: for (i=0; i<count; i++) {
87: SumNorm2(&xin[i*2],&xin[i*2+1],&xout[i*2],&xout[i*2+1]);
88: }
89: } else if (*datatype == MPIU_NORM1_AND_2) {
90: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
91: for (i=0; i<count; i++) {
92: xout[i*3]+= xin[i*3];
93: SumNorm2(&xin[i*3+1],&xin[i*3+2],&xout[i*3+1],&xout[i*3+2]);
94: }
95: } else {
96: (*PetscErrorPrintf)("Can only handle MPIU_NORM* data types");
97: MPI_Abort(MPI_COMM_WORLD,1);
98: }
99: PetscFunctionReturnVoid();
100: }
104: static PetscErrorCode VecNormCompEnd(void)
105: {
109: MPI_Type_free(&MPIU_NORM2);
110: MPI_Type_free(&MPIU_NORM1_AND_2);
111: MPI_Op_free(&MPIU_NORM2_SUM);
112: VecCompInitialized = PETSC_FALSE;
113: return(0);
114: }
118: static PetscErrorCode VecNormCompInit()
119: {
123: MPI_Type_contiguous(sizeof(PetscReal)*2,MPI_BYTE,&MPIU_NORM2);
124: MPI_Type_commit(&MPIU_NORM2);
125: MPI_Type_contiguous(sizeof(PetscReal)*3,MPI_BYTE,&MPIU_NORM1_AND_2);
126: MPI_Type_commit(&MPIU_NORM1_AND_2);
127: MPI_Op_create(SlepcSumNorm2_Local,1,&MPIU_NORM2_SUM);
128: PetscRegisterFinalize(VecNormCompEnd);
129: return(0);
130: }
131: #endif
135: PetscErrorCode VecDestroy_Comp(Vec v)
136: {
137: Vec_Comp *vs = (Vec_Comp*)v->data;
138: PetscInt i;
142: #if defined(PETSC_USE_LOG)
143: PetscLogObjectState((PetscObject)v,"Length=%D",v->map->n);
144: #endif
145: for (i=0;i<vs->nx;i++) {
146: VecDestroy(&vs->x[i]);
147: }
148: if (--vs->n->friends <= 0) {
149: PetscFree(vs->n);
150: }
151: PetscFree(vs->x);
152: PetscFree(vs);
153: return(0);
154: }
156: static struct _VecOps DvOps = {VecDuplicate_Comp, /* 1 */
157: VecDuplicateVecs_Comp,
158: VecDestroyVecs_Comp,
159: VecDot_Comp_MPI,
160: VecMDot_Comp_MPI,
161: VecNorm_Comp_MPI,
162: VecTDot_Comp_MPI,
163: VecMTDot_Comp_MPI,
164: VecScale_Comp,
165: VecCopy_Comp, /* 10 */
166: VecSet_Comp,
167: VecSwap_Comp,
168: VecAXPY_Comp,
169: VecAXPBY_Comp,
170: VecMAXPY_Comp,
171: VecAYPX_Comp,
172: VecWAXPY_Comp,
173: VecAXPBYPCZ_Comp,
174: VecPointwiseMult_Comp,
175: VecPointwiseDivide_Comp,
176: 0, /* 20 */
177: 0,0,
178: 0 /*VecGetArray_Seq*/,
179: VecGetSize_Comp,
180: VecGetLocalSize_Comp,
181: 0/*VecRestoreArray_Seq*/,
182: VecMax_Comp,
183: VecMin_Comp,
184: VecSetRandom_Comp,
185: 0, /* 30 */
186: 0,
187: VecDestroy_Comp,
188: VecView_Comp,
189: 0/*VecPlaceArray_Seq*/,
190: 0/*VecReplaceArray_Seq*/,
191: VecDot_Comp_Seq,
192: 0,
193: VecNorm_Comp_Seq,
194: VecMDot_Comp_Seq,
195: 0, /* 40 */
196: 0,
197: VecReciprocal_Comp,
198: VecConjugate_Comp,
199: 0,0,
200: 0/*VecResetArray_Seq*/,
201: 0,
202: VecMaxPointwiseDivide_Comp,
203: VecPointwiseMax_Comp,
204: VecPointwiseMaxAbs_Comp,
205: VecPointwiseMin_Comp,
206: 0,
207: VecSqrtAbs_Comp,
208: VecAbs_Comp,
209: VecExp_Comp,
210: VecLog_Comp,
211: 0/*VecShift_Comp*/,
212: 0,
213: 0,
214: 0,
215: VecDotNorm2_Comp_MPI
216: };
220: PetscErrorCode VecDuplicateVecs_Comp(Vec w,PetscInt m,Vec *V[])
221: {
223: PetscInt i;
228: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
229: PetscMalloc(m*sizeof(Vec*),V);
230: for (i=0;i<m;i++) { VecDuplicate(w,*V+i); }
231: return(0);
232: }
236: PetscErrorCode VecDestroyVecs_Comp(PetscInt m,Vec v[])
237: {
239: PetscInt i;
243: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
244: for (i=0;i<m;i++) { VecDestroy(&v[i]); }
245: PetscFree(v);
246: return(0);
247: }
251: static PetscErrorCode VecCreate_Comp_Private(Vec v,Vec *x,PetscInt nx,PetscBool x_to_me,Vec_Comp_N *n)
252: {
253: Vec_Comp *s;
255: PetscInt N=0,lN=0,i,k;
258: #if !defined(PETSC_USE_DYNAMIC_LIBRARIES)
259: if (!VecCompInitialized) {
260: VecCompInitialized = PETSC_TRUE;
261: VecRegister(VECCOMP,VecCreate_Comp);
262: VecNormCompInit();
263: }
264: #endif
265: /* Allocate a new Vec_Comp */
266: if (v->data) { PetscFree(v->data); }
267: PetscNewLog(v,Vec_Comp,&s);
268: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
269: v->data = (void*)s;
270: v->petscnative = PETSC_FALSE;
272: /* Allocate the array of Vec, if it is needed to be done */
273: if (x_to_me != PETSC_TRUE) {
274: PetscMalloc(sizeof(Vec)*nx,&s->x);
275: PetscMemcpy(s->x,x,sizeof(Vec)*nx);
276: } else s->x = x;
278: s->nx = nx;
279: for (i=0;i<nx;i++) {
280: VecGetSize(x[i],&k);
281: N+= k;
282: VecGetLocalSize(x[i],&k);
283: lN+= k;
284: }
286: /* Allocate the shared structure, if it is not given */
287: if (!n) {
288: PetscNewLog(v,Vec_Comp_N,&n);
289: s->n = n;
290: n->n = nx;
291: n->N = N;
292: n->lN = lN;
293: n->friends = 1;
294: } else { /* If not, check in the vector in the shared structure */
295: s->n = n;
296: s->n->friends++;
297: s->n->n = nx;
298: }
300: /* Set the virtual sizes as the real sizes of the vector */
301: VecSetSizes(v,s->n->lN,s->n->N);
303: PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
304: return(0);
305: }
309: PETSC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
310: {
314: VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
315: return(0);
316: }
320: /*@C
321: VecCreateComp - Creates a new vector containing several subvectors, each stored separately
323: Collective on Vec
325: Input Parameter:
326: + comm - communicator for the new Vec
327: . Nx - array of (initial) global sizes of child vectors
328: . n - number of child vectors
329: . t - type of the child vectors
330: - Vparent - (optional) template vector
332: Output Parameter:
333: . V - new vector
335: Notes:
336: This is similar to PETSc's VecNest but customized for SLEPc's needs. In particular,
337: the number of child vectors can be modified dynamically, with VecCompSetSubVecs().
339: Level: developer
341: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
342: @*/
343: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
344: {
346: Vec *x;
347: PetscInt i;
350: VecCreate(comm,V);
351: PetscMalloc(n*sizeof(Vec),&x);
352: PetscLogObjectMemory(*V,n*sizeof(Vec));
353: for (i=0;i<n;i++) {
354: VecCreate(comm,&x[i]);
355: VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
356: VecSetType(x[i],t);
357: }
358: VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,
359: Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
360: return(0);
361: }
365: /*@C
366: VecCreateCompWithVecs - Creates a new vector containing several subvectors,
367: each stored separately, from an array of Vecs
369: Collective on Vec
371: Input Parameter:
372: + x - array of Vecs
373: . n - number of child vectors
374: - Vparent - (optional) template vector
376: Output Parameter:
377: . V - new vector
379: Level: developer
381: .seealso: VecCreateComp()
382: @*/
383: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
384: {
386: PetscInt i;
389: VecCreate(PetscObjectComm((PetscObject)x[0]),V);
390: for (i=0;i<n;i++) {
391: PetscObjectReference((PetscObject)x[i]);
392: }
393: VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,
394: Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
395: return(0);
396: }
400: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
401: {
403: Vec *x;
404: PetscInt i;
405: Vec_Comp *s = (Vec_Comp*)win->data;
408: SlepcValidVecComp(win);
409: VecCreate(PetscObjectComm((PetscObject)win),V);
410: PetscMalloc(s->nx*sizeof(Vec),&x);
411: PetscLogObjectMemory(*V,s->nx*sizeof(Vec));
412: for (i=0;i<s->nx;i++) {
413: VecDuplicate(s->x[i],&x[i]);
414: }
415: VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
416: return(0);
417: }
421: /*@C
422: VecCompGetSubVecs - Returns the entire array of vectors defining a compound vector
424: Collective on Vec
426: Input Parameter:
427: . win - compound vector
429: Output Parameter:
430: + n - number of child vectors
431: - x - array of child vectors
433: Level: developer
435: .seealso: VecCreateComp()
436: @*/
437: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
438: {
439: Vec_Comp *s = (Vec_Comp*)win->data;
442: SlepcValidVecComp(win);
443: if (x) *x = s->x;
444: if (n) *n = s->nx;
445: return(0);
446: }
450: /*@C
451: VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
452: of replaces the subvectors
454: Collective on Vec
456: Input Parameters:
457: + win - compound vector
458: . n - number of child vectors
459: - x - array of child vectors
461: Level: developer
463: .seealso: VecCreateComp()
464: @*/
465: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
466: {
467: Vec_Comp *s = (Vec_Comp*)win->data;
471: SlepcValidVecComp(win);
472: if (x) {
473: if (n > s->nx) {
474: PetscFree(s->x);
475: PetscMalloc(sizeof(Vec)*n,&s->x);
476: }
477: PetscMemcpy(s->x,x,sizeof(Vec)*n);
478: s->nx = n;
479: }
480: s->n->n = n;
481: return(0);
482: }
486: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
487: {
489: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
490: PetscInt i;
493: SlepcValidVecComp(v);
494: SlepcValidVecComp(w);
495: for (i=0;i<vs->n->n;i++) {
496: VecAXPY(vs->x[i],alpha,ws->x[i]);
497: }
498: return(0);
499: }
503: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
504: {
506: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
507: PetscInt i;
510: SlepcValidVecComp(v);
511: SlepcValidVecComp(w);
512: for (i=0;i<vs->n->n;i++) {
513: VecAYPX(vs->x[i],alpha,ws->x[i]);
514: }
515: return(0);
516: }
520: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
521: {
523: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
524: PetscInt i;
527: SlepcValidVecComp(v);
528: SlepcValidVecComp(w);
529: for (i=0;i<vs->n->n;i++) {
530: VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
531: }
532: return(0);
533: }
537: PetscErrorCode VecMAXPY_Comp(Vec v,PetscInt n,const PetscScalar *alpha,Vec *w)
538: {
540: Vec_Comp *vs = (Vec_Comp*)v->data;
541: Vec *wx;
542: PetscInt i,j;
545: SlepcValidVecComp(v);
546: for (i=0;i<n;i++) SlepcValidVecComp(w[i]);
548: PetscMalloc(sizeof(Vec)*n,&wx);
550: for (j=0;j<vs->n->n;j++) {
551: for (i=0;i<n;i++) wx[i] = ((Vec_Comp*)w[i]->data)->x[j];
552: VecMAXPY(vs->x[j],n,alpha,wx);
553: }
555: PetscFree(wx);
556: return(0);
557: }
561: PetscErrorCode VecWAXPY_Comp(Vec v,PetscScalar alpha,Vec w,Vec z)
562: {
564: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
565: PetscInt i;
568: SlepcValidVecComp(v);
569: SlepcValidVecComp(w);
570: SlepcValidVecComp(z);
571: for (i=0;i<vs->n->n;i++) {
572: VecWAXPY(vs->x[i],alpha,ws->x[i],zs->x[i]);
573: }
574: return(0);
575: }
579: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
580: {
581: PetscErrorCode ierr;
582: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
583: PetscInt i;
586: SlepcValidVecComp(v);
587: SlepcValidVecComp(w);
588: SlepcValidVecComp(z);
589: for (i=0;i<vs->n->n;i++) {
590: VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
591: }
592: return(0);
593: }
597: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
598: {
599: Vec_Comp *vs = (Vec_Comp*)v->data;
602: SlepcValidVecComp(v);
604: *size = vs->n->N;
605: return(0);
606: }
610: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
611: {
612: Vec_Comp *vs = (Vec_Comp*)v->data;
615: SlepcValidVecComp(v);
617: *size = vs->n->lN;
618: return(0);
619: }
623: PetscErrorCode VecMax_Comp(Vec v,PetscInt *idx,PetscReal *z)
624: {
626: Vec_Comp *vs = (Vec_Comp*)v->data;
627: PetscInt idxp,s=0,s0;
628: PetscReal zp,z0;
629: PetscInt i;
632: SlepcValidVecComp(v);
633: if (!idx && !z) return(0);
635: if (vs->n->n > 0) {
636: VecMax(vs->x[0],idx?&idxp:NULL,&zp);
637: }
638: for (i=1;i<vs->n->n;i++) {
639: VecGetSize(vs->x[i-1],&s0);
640: s+= s0;
641: VecMax(vs->x[i],idx?&idxp:NULL,&z0);
642: if (zp < z0) {
643: if (idx) *idx = s+idxp;
644: zp = z0;
645: }
646: }
647: if (z) *z = zp;
648: return(0);
649: }
653: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
654: {
656: Vec_Comp *vs = (Vec_Comp*)v->data;
657: PetscInt idxp,s=0,s0;
658: PetscReal zp,z0;
659: PetscInt i;
662: SlepcValidVecComp(v);
663: if (!idx && !z) return(0);
665: if (vs->n->n > 0) {
666: VecMin(vs->x[0],idx?&idxp:NULL,&zp);
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)