Actual source code: veccomp.c
slepc-3.8.0 2017-10-20
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2017, 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 MPI_Op MPIU_NORM2_SUM=0;
16: static PetscBool VecCompInitialized = PETSC_FALSE;
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: static void SlepcSumNorm2_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
65: {
66: PetscInt i,count = *cnt;
69: if (*datatype == MPIU_NORM2) {
70: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
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: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
76: for (i=0; i<count; i++) {
77: xout[i*3]+= xin[i*3];
78: SumNorm2(&xin[i*3+1],&xin[i*3+2],&xout[i*3+1],&xout[i*3+2]);
79: }
80: } else {
81: (*PetscErrorPrintf)("Can only handle MPIU_NORM* data types");
82: MPI_Abort(MPI_COMM_WORLD,1);
83: }
84: PetscFunctionReturnVoid();
85: }
87: static PetscErrorCode VecNormCompEnd(void)
88: {
92: MPI_Type_free(&MPIU_NORM2);
93: MPI_Type_free(&MPIU_NORM1_AND_2);
94: MPI_Op_free(&MPIU_NORM2_SUM);
95: VecCompInitialized = PETSC_FALSE;
96: return(0);
97: }
99: static PetscErrorCode VecNormCompInit()
100: {
104: MPI_Type_contiguous(sizeof(PetscReal)*2,MPI_BYTE,&MPIU_NORM2);
105: MPI_Type_commit(&MPIU_NORM2);
106: MPI_Type_contiguous(sizeof(PetscReal)*3,MPI_BYTE,&MPIU_NORM1_AND_2);
107: MPI_Type_commit(&MPIU_NORM1_AND_2);
108: MPI_Op_create(SlepcSumNorm2_Local,1,&MPIU_NORM2_SUM);
109: PetscRegisterFinalize(VecNormCompEnd);
110: return(0);
111: }
113: PetscErrorCode VecDestroy_Comp(Vec v)
114: {
115: Vec_Comp *vs = (Vec_Comp*)v->data;
116: PetscInt i;
120: #if defined(PETSC_USE_LOG)
121: PetscLogObjectState((PetscObject)v,"Length=%D",v->map->n);
122: #endif
123: for (i=0;i<vs->nx;i++) {
124: VecDestroy(&vs->x[i]);
125: }
126: if (--vs->n->friends <= 0) {
127: PetscFree(vs->n);
128: }
129: PetscFree(vs->x);
130: PetscFree(vs);
131: PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",NULL);
132: PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",NULL);
133: return(0);
134: }
136: static struct _VecOps DvOps = {VecDuplicate_Comp, /* 1 */
137: VecDuplicateVecs_Comp,
138: VecDestroyVecs_Comp,
139: VecDot_Comp_MPI,
140: VecMDot_Comp_MPI,
141: VecNorm_Comp_MPI,
142: VecTDot_Comp_MPI,
143: VecMTDot_Comp_MPI,
144: VecScale_Comp,
145: VecCopy_Comp, /* 10 */
146: VecSet_Comp,
147: VecSwap_Comp,
148: VecAXPY_Comp,
149: VecAXPBY_Comp,
150: VecMAXPY_Comp,
151: VecAYPX_Comp,
152: VecWAXPY_Comp,
153: VecAXPBYPCZ_Comp,
154: VecPointwiseMult_Comp,
155: VecPointwiseDivide_Comp,
156: 0, /* 20 */
157: 0,0,
158: 0 /*VecGetArray_Seq*/,
159: VecGetSize_Comp,
160: VecGetLocalSize_Comp,
161: 0/*VecRestoreArray_Seq*/,
162: VecMax_Comp,
163: VecMin_Comp,
164: VecSetRandom_Comp,
165: 0, /* 30 */
166: 0,
167: VecDestroy_Comp,
168: VecView_Comp,
169: 0/*VecPlaceArray_Seq*/,
170: 0/*VecReplaceArray_Seq*/,
171: VecDot_Comp_Seq,
172: VecTDot_Comp_Seq,
173: VecNorm_Comp_Seq,
174: VecMDot_Comp_Seq,
175: VecMTDot_Comp_Seq, /* 40 */
176: 0,
177: VecReciprocal_Comp,
178: VecConjugate_Comp,
179: 0,0,
180: 0/*VecResetArray_Seq*/,
181: 0,
182: VecMaxPointwiseDivide_Comp,
183: VecPointwiseMax_Comp,
184: VecPointwiseMaxAbs_Comp,
185: VecPointwiseMin_Comp,
186: 0,
187: VecSqrtAbs_Comp,
188: VecAbs_Comp,
189: VecExp_Comp,
190: VecLog_Comp,
191: VecShift_Comp,
192: 0,
193: 0,
194: 0,
195: VecDotNorm2_Comp_MPI
196: };
198: PetscErrorCode VecDuplicateVecs_Comp(Vec w,PetscInt m,Vec *V[])
199: {
201: PetscInt i;
206: if (m<=0) SETERRQ1(PetscObjectComm((PetscObject)w),PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
207: PetscMalloc1(m,V);
208: for (i=0;i<m;i++) { VecDuplicate(w,*V+i); }
209: return(0);
210: }
212: PetscErrorCode VecDestroyVecs_Comp(PetscInt m,Vec v[])
213: {
215: PetscInt i;
219: if (m<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
220: for (i=0;i<m;i++) { VecDestroy(&v[i]); }
221: PetscFree(v);
222: return(0);
223: }
225: static PetscErrorCode VecCreate_Comp_Private(Vec v,Vec *x,PetscInt nx,PetscBool x_to_me,Vec_Comp_N *n)
226: {
227: Vec_Comp *s;
229: PetscInt N=0,lN=0,i,k;
232: if (!VecCompInitialized) {
233: VecCompInitialized = PETSC_TRUE;
234: VecRegister(VECCOMP,VecCreate_Comp);
235: VecNormCompInit();
236: }
238: /* Allocate a new Vec_Comp */
239: if (v->data) { PetscFree(v->data); }
240: PetscNewLog(v,&s);
241: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
242: v->data = (void*)s;
243: v->petscnative = PETSC_FALSE;
245: /* Allocate the array of Vec, if it is needed to be done */
246: if (!x_to_me) {
247: if (nx) { PetscMalloc1(nx,&s->x); }
248: if (x) { PetscMemcpy(s->x,x,sizeof(Vec)*nx); }
249: } else s->x = x;
251: s->nx = nx;
253: if (nx) {
254: /* Allocate the shared structure, if it is not given */
255: if (!n) {
256: for (i=0;i<nx;i++) {
257: VecGetSize(x[i],&k);
258: N+= k;
259: VecGetLocalSize(x[i],&k);
260: lN+= k;
261: }
262: PetscNewLog(v,&n);
263: s->n = n;
264: n->n = nx;
265: n->N = N;
266: n->lN = lN;
267: n->friends = 1;
268: } else { /* If not, check in the vector in the shared structure */
269: s->n = n;
270: s->n->friends++;
271: }
273: /* Set the virtual sizes as the real sizes of the vector */
274: VecSetSizes(v,s->n->lN,s->n->N);
275: }
277: PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
278: PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",VecCompSetSubVecs_Comp);
279: PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",VecCompGetSubVecs_Comp);
280: return(0);
281: }
283: PETSC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
284: {
288: VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
289: return(0);
290: }
292: /*@C
293: VecCreateComp - Creates a new vector containing several subvectors,
294: each stored separately.
296: Collective on Vec
298: Input Parameters:
299: + comm - communicator for the new Vec
300: . Nx - array of (initial) global sizes of child vectors
301: . n - number of child vectors
302: . t - type of the child vectors
303: - Vparent - (optional) template vector
305: Output Parameter:
306: . V - new vector
308: Notes:
309: This is similar to PETSc's VecNest but customized for SLEPc's needs. In particular,
310: the number of child vectors can be modified dynamically, with VecCompSetSubVecs().
312: Level: developer
314: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
315: @*/
316: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
317: {
319: Vec *x;
320: PetscInt i;
323: VecCreate(comm,V);
324: PetscMalloc1(n,&x);
325: PetscLogObjectMemory((PetscObject)*V,n*sizeof(Vec));
326: for (i=0;i<n;i++) {
327: VecCreate(comm,&x[i]);
328: VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
329: VecSetType(x[i],t);
330: }
331: VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
332: return(0);
333: }
335: /*@C
336: VecCreateCompWithVecs - Creates a new vector containing several subvectors,
337: each stored separately, from an array of Vecs.
339: Collective on Vec
341: Input Parameters:
342: + x - array of Vecs
343: . n - number of child vectors
344: - Vparent - (optional) template vector
346: Output Parameter:
347: . V - new vector
349: Level: developer
351: .seealso: VecCreateComp()
352: @*/
353: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
354: {
356: PetscInt i;
359: VecCreate(PetscObjectComm((PetscObject)x[0]),V);
360: for (i=0;i<n;i++) {
361: PetscObjectReference((PetscObject)x[i]);
362: }
363: VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
364: return(0);
365: }
367: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
368: {
370: Vec *x;
371: PetscInt i;
372: Vec_Comp *s = (Vec_Comp*)win->data;
375: SlepcValidVecComp(win,1);
376: VecCreate(PetscObjectComm((PetscObject)win),V);
377: PetscMalloc1(s->nx,&x);
378: PetscLogObjectMemory((PetscObject)*V,s->nx*sizeof(Vec));
379: for (i=0;i<s->nx;i++) {
380: if (s->x[i]) {
381: VecDuplicate(s->x[i],&x[i]);
382: } else x[i] = NULL;
383: }
384: VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
385: return(0);
386: }
388: static PetscErrorCode VecCompGetSubVecs_Comp(Vec win,PetscInt *n,const Vec **x)
389: {
390: Vec_Comp *s = (Vec_Comp*)win->data;
393: if (x) *x = s->x;
394: if (n) *n = s->n->n;
395: return(0);
396: }
398: /*@C
399: VecCompGetSubVecs - Returns the entire array of vectors defining a
400: compound vector.
402: Collective on Vec
404: Input Parameter:
405: . win - compound vector
407: Output Parameters:
408: + n - number of child vectors
409: - x - array of child vectors
411: Level: developer
413: .seealso: VecCreateComp()
414: @*/
415: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
416: {
421: PetscUseMethod(win,"VecCompGetSubVecs_C",(Vec,PetscInt*,const Vec**),(win,n,x));
422: return(0);
423: }
425: static PetscErrorCode VecCompSetSubVecs_Comp(Vec win,PetscInt n,Vec *x)
426: {
427: Vec_Comp *s = (Vec_Comp*)win->data;
428: PetscInt i,N,nlocal;
429: Vec_Comp_N *nn;
433: if (!s) SETERRQ(PetscObjectComm((PetscObject)win),PETSC_ERR_ORDER,"Must call VecSetSizes first");
434: if (!s->nx) {
435: /* vector has been created via VecCreate+VecSetType+VecSetSizes, so allocate data structures */
436: PetscMalloc1(n,&s->x);
437: PetscLogObjectMemory((PetscObject)win,n*sizeof(Vec));
438: VecGetSize(win,&N);
439: if (N%n) SETERRQ2(PetscObjectComm((PetscObject)win),1,"Global dimension %D is not divisible by %D",N,n);
440: VecGetLocalSize(win,&nlocal);
441: if (nlocal%n) SETERRQ2(PetscObjectComm((PetscObject)win),1,"Local dimension %D is not divisible by %D",nlocal,n);
442: s->nx = n;
443: for (i=0;i<n;i++) {
444: VecCreate(PetscObjectComm((PetscObject)win),&s->x[i]);
445: VecSetSizes(s->x[i],nlocal/n,N/n);
446: VecSetFromOptions(s->x[i]);
447: }
448: if (!s->n) {
449: PetscNewLog(win,&nn);
450: s->n = nn;
451: nn->N = N;
452: nn->lN = nlocal;
453: nn->friends = 1;
454: }
455: } else if (n > s->nx) SETERRQ1(PetscObjectComm((PetscObject)win),PETSC_ERR_SUP,"Number of child vectors cannot be larger than %D",s->nx);
456: if (x) {
457: PetscMemcpy(s->x,x,sizeof(Vec)*n);
458: }
459: s->n->n = n;
460: return(0);
461: }
463: /*@C
464: VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
465: or replaces the subvectors.
467: Collective on Vec
469: Input Parameters:
470: + win - compound vector
471: . n - number of child vectors
472: - x - array of child vectors
474: Note:
475: It is not possible to increase the number of subvectors with respect to the
476: number set at its creation.
478: Level: developer
480: .seealso: VecCreateComp(), VecCompGetSubVecs()
481: @*/
482: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
483: {
489: PetscTryMethod(win,"VecCompSetSubVecs_C",(Vec,PetscInt,Vec*),(win,n,x));
490: return(0);
491: }
493: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
494: {
496: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
497: PetscInt i;
500: SlepcValidVecComp(v,1);
501: SlepcValidVecComp(w,3);
502: for (i=0;i<vs->n->n;i++) {
503: VecAXPY(vs->x[i],alpha,ws->x[i]);
504: }
505: return(0);
506: }
508: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
509: {
511: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
512: PetscInt i;
515: SlepcValidVecComp(v,1);
516: SlepcValidVecComp(w,3);
517: for (i=0;i<vs->n->n;i++) {
518: VecAYPX(vs->x[i],alpha,ws->x[i]);
519: }
520: return(0);
521: }
523: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
524: {
526: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
527: PetscInt i;
530: SlepcValidVecComp(v,1);
531: SlepcValidVecComp(w,4);
532: for (i=0;i<vs->n->n;i++) {
533: VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
534: }
535: return(0);
536: }
538: PetscErrorCode VecMAXPY_Comp(Vec v,PetscInt n,const PetscScalar *alpha,Vec *w)
539: {
541: Vec_Comp *vs = (Vec_Comp*)v->data;
542: Vec *wx;
543: PetscInt i,j;
546: SlepcValidVecComp(v,1);
547: for (i=0;i<n;i++) SlepcValidVecComp(w[i],4);
549: PetscMalloc1(n,&wx);
551: for (j=0;j<vs->n->n;j++) {
552: for (i=0;i<n;i++) wx[i] = ((Vec_Comp*)w[i]->data)->x[j];
553: VecMAXPY(vs->x[j],n,alpha,wx);
554: }
556: PetscFree(wx);
557: return(0);
558: }
560: PetscErrorCode VecWAXPY_Comp(Vec v,PetscScalar alpha,Vec w,Vec z)
561: {
563: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
564: PetscInt i;
567: SlepcValidVecComp(v,1);
568: SlepcValidVecComp(w,3);
569: SlepcValidVecComp(z,4);
570: for (i=0;i<vs->n->n;i++) {
571: VecWAXPY(vs->x[i],alpha,ws->x[i],zs->x[i]);
572: }
573: return(0);
574: }
576: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
577: {
578: PetscErrorCode ierr;
579: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
580: PetscInt i;
583: SlepcValidVecComp(v,1);
584: SlepcValidVecComp(w,5);
585: SlepcValidVecComp(z,6);
586: for (i=0;i<vs->n->n;i++) {
587: VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
588: }
589: return(0);
590: }
592: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
593: {
594: Vec_Comp *vs = (Vec_Comp*)v->data;
598: if (vs->n) {
599: SlepcValidVecComp(v,1);
600: *size = vs->n->N;
601: } else *size = v->map->N;
602: return(0);
603: }
605: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
606: {
607: Vec_Comp *vs = (Vec_Comp*)v->data;
611: if (vs->n) {
612: SlepcValidVecComp(v,1);
613: *size = vs->n->lN;
614: } else *size = v->map->n;
615: return(0);
616: }
618: PetscErrorCode VecMax_Comp(Vec v,PetscInt *idx,PetscReal *z)
619: {
621: Vec_Comp *vs = (Vec_Comp*)v->data;
622: PetscInt idxp,s=0,s0;
623: PetscReal zp,z0;
624: PetscInt i;
627: SlepcValidVecComp(v,1);
628: if (!idx && !z) return(0);
630: if (vs->n->n > 0) {
631: VecMax(vs->x[0],idx?&idxp:NULL,&zp);
632: } else {
633: zp = PETSC_MIN_REAL;
634: if (idx) idxp = -1;
635: }
636: for (i=1;i<vs->n->n;i++) {
637: VecGetSize(vs->x[i-1],&s0);
638: s += s0;
639: VecMax(vs->x[i],idx?&idxp:NULL,&z0);
640: if (zp < z0) {
641: if (idx) *idx = s+idxp;
642: zp = z0;
643: }
644: }
645: if (z) *z = zp;
646: return(0);
647: }
649: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
650: {
652: Vec_Comp *vs = (Vec_Comp*)v->data;
653: PetscInt idxp,s=0,s0;
654: PetscReal zp,z0;
655: PetscInt i;
658: SlepcValidVecComp(v,1);
659: if (!idx && !z) return(0);
661: if (vs->n->n > 0) {
662: VecMin(vs->x[0],idx?&idxp:NULL,&zp);
663: } else {
664: zp = PETSC_MAX_REAL;
665: if (idx) idxp = -1;
666: }
667: for (i=1;i<vs->n->n;i++) {
668: VecGetSize(vs->x[i-1],&s0);
669: s += s0;
670: VecMin(vs->x[i],idx?&idxp:NULL,&z0);
671: if (zp > z0) {
672: if (idx) *idx = s+idxp;
673: zp = z0;
674: }
675: }
676: if (z) *z = zp;
677: return(0);
678: }
680: PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v,Vec w,PetscReal *m)
681: {
683: Vec_Comp *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
684: PetscReal work;
685: PetscInt i;
688: SlepcValidVecComp(v,1);
689: SlepcValidVecComp(w,2);
690: if (!m || vs->n->n == 0) return(0);
691: VecMaxPointwiseDivide(vs->x[0],ws->x[0],m);
692: for (i=1;i<vs->n->n;i++) {
693: VecMaxPointwiseDivide(vs->x[i],ws->x[i],&work);
694: *m = PetscMax(*m,work);
695: }
696: return(0);
697: }
705: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v) \
706: { \
707: PetscErrorCode ierr; \
708: Vec_Comp *vs = (Vec_Comp*)v->data; \
709: PetscInt i; \
710: \
712: SlepcValidVecComp(v,1); \
713: for (i=0;i<vs->n->n;i++) { \
714: __COMPOSE2__(Vec,NAME)(vs->x[i]); \
715: } \
716: return(0);\
717: }
719: __FUNC_TEMPLATE1__(Conjugate)
720: __FUNC_TEMPLATE1__(Reciprocal)
721: __FUNC_TEMPLATE1__(SqrtAbs)
722: __FUNC_TEMPLATE1__(Abs)
723: __FUNC_TEMPLATE1__(Exp)
724: __FUNC_TEMPLATE1__(Log)
727: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,T0 __a) \
728: { \
729: PetscErrorCode ierr; \
730: Vec_Comp *vs = (Vec_Comp*)v->data; \
731: PetscInt i; \
732: \
734: SlepcValidVecComp(v,1); \
735: for (i=0;i<vs->n->n;i++) { \
736: __COMPOSE2__(Vec,NAME)(vs->x[i],__a); \
737: } \
738: return(0);\
739: }
741: __FUNC_TEMPLATE2__(Set,PetscScalar)
742: __FUNC_TEMPLATE2__(View,PetscViewer)
743: __FUNC_TEMPLATE2__(Scale,PetscScalar)
744: __FUNC_TEMPLATE2__(SetRandom,PetscRandom)
745: __FUNC_TEMPLATE2__(Shift,PetscScalar)
748: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w) \
749: { \
750: PetscErrorCode ierr; \
751: Vec_Comp *vs = (Vec_Comp*)v->data,\
752: *ws = (Vec_Comp*)w->data; \
753: PetscInt i; \
754: \
756: SlepcValidVecComp(v,1); \
757: SlepcValidVecComp(w,2); \
758: for (i=0;i<vs->n->n;i++) { \
759: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i]); \
760: } \
761: return(0);\
762: }
764: __FUNC_TEMPLATE3__(Copy)
765: __FUNC_TEMPLATE3__(Swap)
768: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w,Vec z) \
769: { \
770: PetscErrorCode ierr; \
771: Vec_Comp *vs = (Vec_Comp*)v->data, \
772: *ws = (Vec_Comp*)w->data, \
773: *zs = (Vec_Comp*)z->data; \
774: PetscInt i; \
775: \
777: SlepcValidVecComp(v,1); \
778: SlepcValidVecComp(w,2); \
779: SlepcValidVecComp(z,3); \
780: for (i=0;i<vs->n->n;i++) { \
781: __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i],zs->x[i]); \
782: } \
783: return(0);\
784: }
786: __FUNC_TEMPLATE4__(PointwiseMax)
787: __FUNC_TEMPLATE4__(PointwiseMaxAbs)
788: __FUNC_TEMPLATE4__(PointwiseMin)
789: __FUNC_TEMPLATE4__(PointwiseMult)
790: __FUNC_TEMPLATE4__(PointwiseDivide)