Actual source code: veccomp.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #include <slepc/private/vecimplslepc.h> /*I "slepcvec.h" I*/
13: /* Private MPI datatypes and operators */
14: static MPI_Datatype MPIU_NORM2=0, MPIU_NORM1_AND_2=0;
15: static PetscBool VecCompInitialized = PETSC_FALSE;
16: MPI_Op MPIU_NORM2_SUM=0;
18: /* Private functions */
19: PETSC_STATIC_INLINE void SumNorm2(PetscReal*,PetscReal*,PetscReal*,PetscReal*);
20: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal,PetscReal);
21: PETSC_STATIC_INLINE void AddNorm2(PetscReal*,PetscReal*,PetscReal);
22: static PetscErrorCode VecCompSetSubVecs_Comp(Vec,PetscInt,Vec*);
23: static PetscErrorCode VecCompGetSubVecs_Comp(Vec,PetscInt*,const Vec**);
25: #include "veccomp0.h"
28: #include "veccomp0.h"
30: PETSC_STATIC_INLINE void SumNorm2(PetscReal *ssq0,PetscReal *scale0,PetscReal *ssq1,PetscReal *scale1)
31: {
32: PetscReal q;
33: if (*scale0 > *scale1) {
34: q = *scale1/(*scale0);
35: *ssq1 = *ssq0 + q*q*(*ssq1);
36: *scale1 = *scale0;
37: } else {
38: q = *scale0/(*scale1);
39: *ssq1 += q*q*(*ssq0);
40: }
41: }
43: PETSC_STATIC_INLINE PetscReal GetNorm2(PetscReal ssq,PetscReal scale)
44: {
45: return scale*PetscSqrtReal(ssq);
46: }
48: PETSC_STATIC_INLINE void AddNorm2(PetscReal *ssq,PetscReal *scale,PetscReal x)
49: {
50: PetscReal absx,q;
51: if (x != 0.0) {
52: absx = PetscAbs(x);
53: if (*scale < absx) {
54: q = *scale/absx;
55: *ssq = 1.0 + *ssq*q*q;
56: *scale = absx;
57: } else {
58: q = absx/(*scale);
59: *ssq += q*q;
60: }
61: }
62: }
64: SLEPC_EXTERN void SlepcSumNorm2_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
65: {
66: PetscInt i,count = *cnt;
67: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
70: if (*datatype == MPIU_NORM2) {
71: for (i=0;i<count;i++) {
72: SumNorm2(&xin[i*2],&xin[i*2+1],&xout[i*2],&xout[i*2+1]);
73: }
74: } else if (*datatype == MPIU_NORM1_AND_2) {
75: for (i=0;i<count;i++) {
76: xout[i*3] += xin[i*3];
77: SumNorm2(&xin[i*3+1],&xin[i*3+2],&xout[i*3+1],&xout[i*3+2]);
78: }
79: } else {
80: (*PetscErrorPrintf)("Can only handle MPIU_NORM* data types");
81: MPI_Abort(MPI_COMM_WORLD,1);
82: }
83: PetscFunctionReturnVoid();
84: }
86: static PetscErrorCode VecCompNormEnd(void)
87: {
91: MPI_Type_free(&MPIU_NORM2);
92: MPI_Type_free(&MPIU_NORM1_AND_2);
93: MPI_Op_free(&MPIU_NORM2_SUM);
94: VecCompInitialized = PETSC_FALSE;
95: return(0);
96: }
98: static PetscErrorCode VecCompNormInit(void)
99: {
103: MPI_Type_contiguous(2,MPIU_REAL,&MPIU_NORM2);
104: MPI_Type_commit(&MPIU_NORM2);
105: MPI_Type_contiguous(3,MPIU_REAL,&MPIU_NORM1_AND_2);
106: MPI_Type_commit(&MPIU_NORM1_AND_2);
107: MPI_Op_create(SlepcSumNorm2_Local,PETSC_TRUE,&MPIU_NORM2_SUM);
108: PetscRegisterFinalize(VecCompNormEnd);
109: return(0);
110: }
112: PetscErrorCode VecDestroy_Comp(Vec v)
113: {
114: Vec_Comp *vs = (Vec_Comp*)v->data;
115: PetscInt i;
119: #if defined(PETSC_USE_LOG)
120: PetscLogObjectState((PetscObject)v,"Length=%D",v->map->n);
121: #endif
122: for (i=0;i<vs->nx;i++) {
123: VecDestroy(&vs->x[i]);
124: }
125: if (--vs->n->friends <= 0) {
126: PetscFree(vs->n);
127: }
128: PetscFree(vs->x);
129: PetscFree(vs);
130: PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",NULL);
131: PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",NULL);
132: return(0);
133: }
135: static struct _VecOps DvOps = {VecDuplicate_Comp, /* 1 */
136: VecDuplicateVecs_Comp,
137: VecDestroyVecs_Comp,
138: VecDot_Comp_MPI,
139: VecMDot_Comp_MPI,
140: VecNorm_Comp_MPI,
141: VecTDot_Comp_MPI,
142: VecMTDot_Comp_MPI,
143: VecScale_Comp,
144: VecCopy_Comp, /* 10 */
145: VecSet_Comp,
146: VecSwap_Comp,
147: VecAXPY_Comp,
148: VecAXPBY_Comp,
149: VecMAXPY_Comp,
150: VecAYPX_Comp,
151: VecWAXPY_Comp,
152: VecAXPBYPCZ_Comp,
153: VecPointwiseMult_Comp,
154: VecPointwiseDivide_Comp,
155: 0, /* 20 */
156: 0,0,
157: 0 /*VecGetArray_Seq*/,
158: VecGetSize_Comp,
159: VecGetLocalSize_Comp,
160: 0/*VecRestoreArray_Seq*/,
161: VecMax_Comp,
162: VecMin_Comp,
163: VecSetRandom_Comp,
164: 0, /* 30 */
165: 0,
166: VecDestroy_Comp,
167: VecView_Comp,
168: 0/*VecPlaceArray_Seq*/,
169: 0/*VecReplaceArray_Seq*/,
170: VecDot_Comp_Seq,
171: VecTDot_Comp_Seq,
172: VecNorm_Comp_Seq,
173: VecMDot_Comp_Seq,
174: VecMTDot_Comp_Seq, /* 40 */
175: 0,
176: VecReciprocal_Comp,
177: VecConjugate_Comp,
178: 0,0,
179: 0/*VecResetArray_Seq*/,
180: 0,
181: VecMaxPointwiseDivide_Comp,
182: VecPointwiseMax_Comp,
183: VecPointwiseMaxAbs_Comp,
184: VecPointwiseMin_Comp,
185: 0,
186: VecSqrtAbs_Comp,
187: VecAbs_Comp,
188: VecExp_Comp,
189: VecLog_Comp,
190: VecShift_Comp,
191: 0,
192: 0,
193: 0,
194: VecDotNorm2_Comp_MPI
195: };
197: PetscErrorCode VecDuplicateVecs_Comp(Vec w,PetscInt m,Vec *V[])
198: {
200: PetscInt i;
205: if (m<=0) SETERRQ1(PetscObjectComm((PetscObject)w),PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
206: PetscMalloc1(m,V);
207: for (i=0;i<m;i++) { VecDuplicate(w,*V+i); }
208: return(0);
209: }
211: PetscErrorCode VecDestroyVecs_Comp(PetscInt m,Vec v[])
212: {
214: PetscInt i;
218: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
219: for (i=0;i<m;i++) { VecDestroy(&v[i]); }
220: PetscFree(v);
221: return(0);
222: }
224: static PetscErrorCode VecCreate_Comp_Private(Vec v,Vec *x,PetscInt nx,PetscBool x_to_me,Vec_Comp_N *n)
225: {
226: Vec_Comp *s;
228: PetscInt N=0,lN=0,i,k;
231: if (!VecCompInitialized) {
232: VecCompInitialized = PETSC_TRUE;
233: VecRegister(VECCOMP,VecCreate_Comp);
234: VecCompNormInit();
235: }
237: /* Allocate a new Vec_Comp */
238: if (v->data) { PetscFree(v->data); }
239: PetscNewLog(v,&s);
240: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
241: v->data = (void*)s;
242: v->petscnative = PETSC_FALSE;
244: /* Allocate the array of Vec, if it is needed to be done */
245: if (!x_to_me) {
246: if (nx) { PetscMalloc1(nx,&s->x); }
247: if (x) { PetscMemcpy(s->x,x,sizeof(Vec)*nx); }
248: } else s->x = x;
250: s->nx = nx;
252: if (nx) {
253: /* Allocate the shared structure, if it is not given */
254: if (!n) {
255: for (i=0;i<nx;i++) {
256: VecGetSize(x[i],&k);
257: N+= k;
258: VecGetLocalSize(x[i],&k);
259: lN+= k;
260: }
261: PetscNewLog(v,&n);
262: s->n = n;
263: n->n = nx;
264: n->N = N;
265: n->lN = lN;
266: n->friends = 1;
267: } else { /* If not, check in the vector in the shared structure */
268: s->n = n;
269: s->n->friends++;
270: }
272: /* Set the virtual sizes as the real sizes of the vector */
273: VecSetSizes(v,s->n->lN,s->n->N);
274: }
276: PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
277: PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",VecCompSetSubVecs_Comp);
278: PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",VecCompGetSubVecs_Comp);
279: return(0);
280: }
282: SLEPC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
283: {
287: VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
288: return(0);
289: }
291: /*@C
292: VecCreateComp - Creates a new vector containing several subvectors,
293: each stored separately.
295: Collective on Vec
297: Input Parameters:
298: + comm - communicator for the new Vec
299: . Nx - array of (initial) global sizes of child vectors
300: . n - number of child vectors
301: . t - type of the child vectors
302: - Vparent - (optional) template vector
304: Output Parameter:
305: . V - new vector
307: Notes:
308: This is similar to PETSc's VecNest but customized for SLEPc's needs. In particular,
309: the number of child vectors can be modified dynamically, with VecCompSetSubVecs().
311: Level: developer
313: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
314: @*/
315: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
316: {
318: Vec *x;
319: PetscInt i;
322: VecCreate(comm,V);
323: PetscMalloc1(n,&x);
324: PetscLogObjectMemory((PetscObject)*V,n*sizeof(Vec));
325: for (i=0;i<n;i++) {
326: VecCreate(comm,&x[i]);
327: VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
328: VecSetType(x[i],t);
329: }
330: VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
331: return(0);
332: }
334: /*@C
335: VecCreateCompWithVecs - Creates a new vector containing several subvectors,
336: each stored separately, from an array of Vecs.
338: Collective on Vec
340: Input Parameters:
341: + x - array of Vecs
342: . n - number of child vectors
343: - Vparent - (optional) template vector
345: Output Parameter:
346: . V - new vector
348: Level: developer
350: .seealso: VecCreateComp()
351: @*/
352: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
353: {
355: PetscInt i;
358: VecCreate(PetscObjectComm((PetscObject)x[0]),V);
359: for (i=0;i<n;i++) {
360: PetscObjectReference((PetscObject)x[i]);
361: }
362: VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
363: return(0);
364: }
366: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
367: {
369: Vec *x;
370: PetscInt i;
371: Vec_Comp *s = (Vec_Comp*)win->data;
374: SlepcValidVecComp(win,1);
375: VecCreate(PetscObjectComm((PetscObject)win),V);
376: PetscMalloc1(s->nx,&x);
377: PetscLogObjectMemory((PetscObject)*V,s->nx*sizeof(Vec));
378: for (i=0;i<s->nx;i++) {
379: if (s->x[i]) {
380: VecDuplicate(s->x[i],&x[i]);
381: } else x[i] = NULL;
382: }
383: VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
384: return(0);
385: }
387: static PetscErrorCode VecCompGetSubVecs_Comp(Vec win,PetscInt *n,const Vec **x)
388: {
389: Vec_Comp *s = (Vec_Comp*)win->data;
392: if (x) *x = s->x;
393: if (n) *n = s->n->n;
394: return(0);
395: }
397: /*@C
398: VecCompGetSubVecs - Returns the entire array of vectors defining a
399: compound vector.
401: Collective on Vec
403: Input Parameter:
404: . win - compound vector
406: Output Parameters:
407: + n - number of child vectors
408: - x - array of child vectors
410: Level: developer
412: .seealso: VecCreateComp()
413: @*/
414: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
415: {
420: PetscUseMethod(win,"VecCompGetSubVecs_C",(Vec,PetscInt*,const Vec**),(win,n,x));
421: return(0);
422: }
424: static PetscErrorCode VecCompSetSubVecs_Comp(Vec win,PetscInt n,Vec *x)
425: {
426: Vec_Comp *s = (Vec_Comp*)win->data;
427: PetscInt i,N,nlocal;
428: Vec_Comp_N *nn;
432: if (!s) SETERRQ(PetscObjectComm((PetscObject)win),PETSC_ERR_ORDER,"Must call VecSetSizes first");
433: if (!s->nx) {
434: /* vector has been created via VecCreate+VecSetType+VecSetSizes, so allocate data structures */
435: PetscMalloc1(n,&s->x);
436: PetscLogObjectMemory((PetscObject)win,n*sizeof(Vec));
437: VecGetSize(win,&N);
438: if (N%n) SETERRQ2(PetscObjectComm((PetscObject)win),1,"Global dimension %D is not divisible by %D",N,n);
439: VecGetLocalSize(win,&nlocal);
440: if (nlocal%n) SETERRQ2(PetscObjectComm((PetscObject)win),1,"Local dimension %D is not divisible by %D",nlocal,n);
441: s->nx = n;
442: for (i=0;i<n;i++) {
443: VecCreate(PetscObjectComm((PetscObject)win),&s->x[i]);
444: VecSetSizes(s->x[i],nlocal/n,N/n);
445: VecSetFromOptions(s->x[i]);
446: }
447: if (!s->n) {
448: PetscNewLog(win,&nn);
449: s->n = nn;
450: nn->N = N;
451: nn->lN = nlocal;
452: nn->friends = 1;
453: }
454: } else if (n > s->nx) SETERRQ1(PetscObjectComm((PetscObject)win),PETSC_ERR_SUP,"Number of child vectors cannot be larger than %D",s->nx);
455: if (x) {
456: PetscMemcpy(s->x,x,sizeof(Vec)*n);
457: }
458: s->n->n = n;
459: return(0);
460: }
462: /*@C
463: VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
464: or replaces the subvectors.
466: Collective on Vec
468: Input Parameters:
469: + win - compound vector
470: . n - number of child vectors
471: - x - array of child vectors
473: Note:
474: It is not possible to increase the number of subvectors with respect to the
475: number set at its creation.
477: Level: developer
479: .seealso: VecCreateComp(), VecCompGetSubVecs()
480: @*/
481: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
482: {
488: PetscTryMethod(win,"VecCompSetSubVecs_C",(Vec,PetscInt,Vec*),(win,n,x));
489: return(0);
490: }
492: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
493: {
495: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
496: PetscInt i;
499: SlepcValidVecComp(v,1);
500: SlepcValidVecComp(w,3);
501: for (i=0;i<vs->n->n;i++) {
502: VecAXPY(vs->x[i],alpha,ws->x[i]);
503: }
504: return(0);
505: }
507: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
508: {
510: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
511: PetscInt i;
514: SlepcValidVecComp(v,1);
515: SlepcValidVecComp(w,3);
516: for (i=0;i<vs->n->n;i++) {
517: VecAYPX(vs->x[i],alpha,ws->x[i]);
518: }
519: return(0);
520: }
522: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
523: {
525: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
526: PetscInt i;
529: SlepcValidVecComp(v,1);
530: SlepcValidVecComp(w,4);
531: for (i=0;i<vs->n->n;i++) {
532: VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
533: }
534: return(0);
535: }
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,1);
546: for (i=0;i<n;i++) SlepcValidVecComp(w[i],4);
548: PetscMalloc1(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: }
559: PetscErrorCode VecWAXPY_Comp(Vec v,PetscScalar alpha,Vec w,Vec z)
560: {
562: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
563: PetscInt i;
566: SlepcValidVecComp(v,1);
567: SlepcValidVecComp(w,3);
568: SlepcValidVecComp(z,4);
569: for (i=0;i<vs->n->n;i++) {
570: VecWAXPY(vs->x[i],alpha,ws->x[i],zs->x[i]);
571: }
572: return(0);
573: }
575: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
576: {
577: PetscErrorCode ierr;
578: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
579: PetscInt i;
582: SlepcValidVecComp(v,1);
583: SlepcValidVecComp(w,5);
584: SlepcValidVecComp(z,6);
585: for (i=0;i<vs->n->n;i++) {
586: VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
587: }
588: return(0);
589: }
591: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
592: {
593: Vec_Comp *vs = (Vec_Comp*)v->data;
597: if (vs->n) {
598: SlepcValidVecComp(v,1);
599: *size = vs->n->N;
600: } else *size = v->map->N;
601: return(0);
602: }
604: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
605: {
606: Vec_Comp *vs = (Vec_Comp*)v->data;
610: if (vs->n) {
611: SlepcValidVecComp(v,1);
612: *size = vs->n->lN;
613: } else *size = v->map->n;
614: return(0);
615: }
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,1);
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: }
648: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
649: {
651: Vec_Comp *vs = (Vec_Comp*)v->data;
652: PetscInt idxp,s=0,s0;
653: PetscReal zp,z0;
654: PetscInt i;
657: SlepcValidVecComp(v,1);
658: if (!idx && !z) return(0);
660: if (vs->n->n > 0) {
661: VecMin(vs->x[0],idx?&idxp:NULL,&zp);
662: } else {
663: zp = PETSC_MAX_REAL;
664: if (idx) idxp = -1;
665: }
666: for (i=1;i<vs->n->n;i++) {
667: VecGetSize(vs->x[i-1],&s0);
668: s += s0;
669: VecMin(vs->x[i],idx?&idxp:NULL,&z0);
670: if (zp > z0) {
671: if (idx) *idx = s+idxp;
672: zp = z0;
673: }
674: }
675: if (z) *z = zp;
676: return(0);
677: }
679: PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v,Vec w,PetscReal *m)
680: {
682: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
683: PetscReal work;
684: PetscInt i;
687: SlepcValidVecComp(v,1);
688: SlepcValidVecComp(w,2);
689: if (!m || vs->n->n == 0) return(0);
690: VecMaxPointwiseDivide(vs->x[0],ws->x[0],m);
691: for (i=1;i<vs->n->n;i++) {
692: VecMaxPointwiseDivide(vs->x[i],ws->x[i],&work);
693: *m = PetscMax(*m,work);
694: }
695: return(0);
696: }
704: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v) \
705: { \
706: PetscErrorCode ierr; \
707: Vec_Comp *vs = (Vec_Comp*)v->data; \
708: PetscInt i; \
709: \
711: SlepcValidVecComp(v,1); \
712: for (i=0;i<vs->n->n;i++) { \
713: __COMPOSE2__(Vec,NAME)(vs->x[i]); \
714: } \
715: return(0);\
716: }
718: __FUNC_TEMPLATE1__(Conjugate)
719: __FUNC_TEMPLATE1__(Reciprocal)
720: __FUNC_TEMPLATE1__(SqrtAbs)
721: __FUNC_TEMPLATE1__(Abs)
722: __FUNC_TEMPLATE1__(Exp)
723: __FUNC_TEMPLATE1__(Log)
726: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,T0 __a) \
727: { \
728: PetscErrorCode ierr; \
729: Vec_Comp *vs = (Vec_Comp*)v->data; \
730: PetscInt i; \
731: \
733: SlepcValidVecComp(v,1); \
734: for (i=0;i<vs->n->n;i++) { \
735: __COMPOSE2__(Vec,NAME)(vs->x[i],__a); \
736: } \
737: return(0);\
738: }
740: __FUNC_TEMPLATE2__(Set,PetscScalar)
741: __FUNC_TEMPLATE2__(View,PetscViewer)
742: __FUNC_TEMPLATE2__(Scale,PetscScalar)
743: __FUNC_TEMPLATE2__(SetRandom,PetscRandom)
744: __FUNC_TEMPLATE2__(Shift,PetscScalar)
747: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w) \
748: { \
749: PetscErrorCode ierr; \
750: Vec_Comp *vs = (Vec_Comp*)v->data,\
751: *ws = (Vec_Comp*)w->data; \
752: PetscInt i; \
753: \
755: SlepcValidVecComp(v,1); \
756: SlepcValidVecComp(w,2); \
757: for (i=0;i<vs->n->n;i++) { \
758: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i]); \
759: } \
760: return(0);\
761: }
763: __FUNC_TEMPLATE3__(Copy)
764: __FUNC_TEMPLATE3__(Swap)
767: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w,Vec z) \
768: { \
769: PetscErrorCode ierr; \
770: Vec_Comp *vs = (Vec_Comp*)v->data, \
771: *ws = (Vec_Comp*)w->data, \
772: *zs = (Vec_Comp*)z->data; \
773: PetscInt i; \
774: \
776: SlepcValidVecComp(v,1); \
777: SlepcValidVecComp(w,2); \
778: SlepcValidVecComp(z,3); \
779: for (i=0;i<vs->n->n;i++) { \
780: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i],zs->x[i]); \
781: } \
782: return(0);\
783: }
785: __FUNC_TEMPLATE4__(PointwiseMax)
786: __FUNC_TEMPLATE4__(PointwiseMaxAbs)
787: __FUNC_TEMPLATE4__(PointwiseMin)
788: __FUNC_TEMPLATE4__(PointwiseMult)
789: __FUNC_TEMPLATE4__(PointwiseDivide)