Actual source code: veccomp.c

slepc-3.5.2 2014-10-10
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2014, 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:     N = n->N;
290:     lN = n->lN;
291:     s->n = n;
292:     s->n->friends++;
293:   }

295:   /* Set the virtual sizes as the real sizes of the vector */
296:   VecSetSizes(v,s->n->lN,s->n->N);

298:   PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
299:   return(0);
300: }

304: PETSC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
305: {

309:   VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
310:   return(0);
311: }

315: /*@C
316:    VecCreateComp - Creates a new vector containing several subvectors, each stored separately

318:    Collective on Vec

320:    Input Parameter:
321: +  comm - communicator for the new Vec
322: .  Nx   - array of (initial) global sizes of child vectors
323: .  n    - number of child vectors
324: .  t    - type of the child vectors
325: -  Vparent - (optional) template vector

327:    Output Parameter:
328: .  V - new vector

330:    Notes:
331:    This is similar to PETSc's VecNest but customized for SLEPc's needs. In particular,
332:    the number of child vectors can be modified dynamically, with VecCompSetSubVecs().

334:    Level: developer

336: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
337: @*/
338: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
339: {
341:   Vec            *x;
342:   PetscInt       i;

345:   VecCreate(comm,V);
346:   PetscMalloc1(n,&x);
347:   PetscLogObjectMemory((PetscObject)*V,n*sizeof(Vec));
348:   for (i=0;i<n;i++) {
349:     VecCreate(comm,&x[i]);
350:     VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
351:     VecSetType(x[i],t);
352:   }
353:   VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,
354:                            Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
355:   return(0);
356: }

360: /*@C
361:    VecCreateCompWithVecs - Creates a new vector containing several subvectors,
362:    each stored separately, from an array of Vecs

364:    Collective on Vec

366:    Input Parameter:
367: +  x - array of Vecs
368: .  n - number of child vectors
369: -  Vparent - (optional) template vector

371:    Output Parameter:
372: .  V - new vector

374:    Level: developer

376: .seealso: VecCreateComp()
377: @*/
378: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
379: {
381:   PetscInt       i;

384:   VecCreate(PetscObjectComm((PetscObject)x[0]),V);
385:   for (i=0;i<n;i++) {
386:     PetscObjectReference((PetscObject)x[i]);
387:   }
388:   VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,
389:                            Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
390:   return(0);
391: }

395: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
396: {
398:   Vec            *x;
399:   PetscInt       i;
400:   Vec_Comp       *s = (Vec_Comp*)win->data;

403:   SlepcValidVecComp(win);
404:   VecCreate(PetscObjectComm((PetscObject)win),V);
405:   PetscMalloc1(s->nx,&x);
406:   PetscLogObjectMemory((PetscObject)*V,s->nx*sizeof(Vec));
407:   for (i=0;i<s->nx;i++) {
408:     if (s->x[i]) {
409:       VecDuplicate(s->x[i],&x[i]);
410:     } else {
411:       x[i] = NULL;
412:     }
413:   }
414:   VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
415:   return(0);
416: }

420: /*@C
421:    VecCompGetSubVecs - Returns the entire array of vectors defining a compound vector

423:    Collective on Vec

425:    Input Parameter:
426: .  win - compound vector

428:    Output Parameter:
429: +  n - number of child vectors
430: -  x - array of child vectors

432:    Level: developer

434: .seealso: VecCreateComp()
435: @*/
436: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
437: {
438:   Vec_Comp *s = (Vec_Comp*)win->data;

441:   SlepcValidVecComp(win);
442:   if (x) *x = s->x;
443:   if (n) *n = s->nx;
444:   return(0);
445: }

449: /*@C
450:    VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
451:    of replaces the subvectors

453:    Collective on Vec

455:    Input Parameters:
456: +  win - compound vector
457: .  n - number of child vectors
458: -  x - array of child vectors

460:    Level: developer

462: .seealso: VecCreateComp()
463: @*/
464: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
465: {
466:   Vec_Comp       *s = (Vec_Comp*)win->data;

470:   if (x) {
471:     if (n > s->nx) {
472:       PetscFree(s->x);
473:       PetscMalloc(sizeof(Vec)*n,&s->x);
474:     }
475:     PetscMemcpy(s->x,x,sizeof(Vec)*n);
476:     s->nx = n;
477:   }
478:   s->n->n = n;
479:   return(0);
480: }

484: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
485: {
487:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
488:   PetscInt       i;

491:   SlepcValidVecComp(v);
492:   SlepcValidVecComp(w);
493:   for (i=0;i<vs->n->n;i++) {
494:     VecAXPY(vs->x[i],alpha,ws->x[i]);
495:   }
496:   return(0);
497: }

501: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
502: {
504:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
505:   PetscInt       i;

508:   SlepcValidVecComp(v);
509:   SlepcValidVecComp(w);
510:   for (i=0;i<vs->n->n;i++) {
511:     VecAYPX(vs->x[i],alpha,ws->x[i]);
512:   }
513:   return(0);
514: }

518: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
519: {
521:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
522:   PetscInt       i;

525:   SlepcValidVecComp(v);
526:   SlepcValidVecComp(w);
527:   for (i=0;i<vs->n->n;i++) {
528:     VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
529:   }
530:   return(0);
531: }

535: PetscErrorCode VecMAXPY_Comp(Vec v,PetscInt n,const PetscScalar *alpha,Vec *w)
536: {
538:   Vec_Comp       *vs = (Vec_Comp*)v->data;
539:   Vec            *wx;
540:   PetscInt       i,j;

543:   SlepcValidVecComp(v);
544:   for (i=0;i<n;i++) SlepcValidVecComp(w[i]);

546:   PetscMalloc(sizeof(Vec)*n,&wx);

548:   for (j=0;j<vs->n->n;j++) {
549:     for (i=0;i<n;i++) wx[i] = ((Vec_Comp*)w[i]->data)->x[j];
550:     VecMAXPY(vs->x[j],n,alpha,wx);
551:   }

553:   PetscFree(wx);
554:   return(0);
555: }

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);
567:   SlepcValidVecComp(w);
568:   SlepcValidVecComp(z);
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: }

577: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
578: {
579:   PetscErrorCode  ierr;
580:   Vec_Comp        *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
581:   PetscInt        i;

584:   SlepcValidVecComp(v);
585:   SlepcValidVecComp(w);
586:   SlepcValidVecComp(z);
587:   for (i=0;i<vs->n->n;i++) {
588:     VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
589:   }
590:   return(0);
591: }

595: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
596: {
597:   Vec_Comp *vs = (Vec_Comp*)v->data;

600:   SlepcValidVecComp(v);
602:   *size = vs->n->N;
603:   return(0);
604: }

608: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
609: {
610:   Vec_Comp *vs = (Vec_Comp*)v->data;

613:   SlepcValidVecComp(v);
615:   *size = vs->n->lN;
616:   return(0);
617: }

621: PetscErrorCode VecMax_Comp(Vec v,PetscInt *idx,PetscReal *z)
622: {
624:   Vec_Comp       *vs = (Vec_Comp*)v->data;
625:   PetscInt       idxp,s=0,s0;
626:   PetscReal      zp,z0;
627:   PetscInt       i;

630:   SlepcValidVecComp(v);
631:   if (!idx && !z) return(0);

633:   if (vs->n->n > 0) {
634:     VecMax(vs->x[0],idx?&idxp:NULL,&zp);
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: }

651: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
652: {
654:   Vec_Comp       *vs = (Vec_Comp*)v->data;
655:   PetscInt       idxp,s=0,s0;
656:   PetscReal      zp,z0;
657:   PetscInt       i;

660:   SlepcValidVecComp(v);
661:   if (!idx && !z) return(0);

663:   if (vs->n->n > 0) {
664:     VecMin(vs->x[0],idx?&idxp:NULL,&zp);
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: }

681: PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v,Vec w,PetscReal *m)
682: {
684:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
685:   PetscReal      work;
686:   PetscInt       i;

689:   SlepcValidVecComp(v);
690:   SlepcValidVecComp(w);
691:   if (!m || vs->n->n == 0) return(0);
692:   VecMaxPointwiseDivide(vs->x[0],ws->x[0],m);
693:   for (i=1;i<vs->n->n;i++) {
694:     VecMaxPointwiseDivide(vs->x[i],ws->x[i],&work);
695:     *m = PetscMax(*m,work);
696:   }
697:   return(0);
698: }



706: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v) \
707: { \
708:   PetscErrorCode  ierr; \
709:   Vec_Comp        *vs = (Vec_Comp*)v->data; \
710:   PetscInt        i; \
711: \
713:   SlepcValidVecComp(v); \
714:   for (i=0;i<vs->n->n;i++) { \
715:     __COMPOSE2__(Vec,NAME)(vs->x[i]); \
716:   } \
717:   return(0);\
718: }

722: __FUNC_TEMPLATE1__(Conjugate)

726: __FUNC_TEMPLATE1__(Reciprocal)

730: __FUNC_TEMPLATE1__(SqrtAbs)

734: __FUNC_TEMPLATE1__(Abs)

738: __FUNC_TEMPLATE1__(Exp)

742: __FUNC_TEMPLATE1__(Log)


746: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,T0 __a) \
747: { \
748:   PetscErrorCode  ierr; \
749:   Vec_Comp        *vs = (Vec_Comp*)v->data; \
750:   PetscInt        i; \
751: \
753:   SlepcValidVecComp(v); \
754:   for (i=0;i<vs->n->n;i++) { \
755:     __COMPOSE2__(Vec,NAME)(vs->x[i],__a); \
756:   } \
757:   return(0);\
758: }

762: __FUNC_TEMPLATE2__(Set,PetscScalar)

766: __FUNC_TEMPLATE2__(View,PetscViewer)

770: __FUNC_TEMPLATE2__(Scale,PetscScalar)

774: __FUNC_TEMPLATE2__(SetRandom,PetscRandom)

778: __FUNC_TEMPLATE2__(Shift,PetscScalar)


782: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w) \
783: { \
784:   PetscErrorCode  ierr; \
785:   Vec_Comp        *vs = (Vec_Comp*)v->data,\
786:                   *ws = (Vec_Comp*)w->data; \
787:   PetscInt        i; \
788: \
790:   SlepcValidVecComp(v); \
791:   SlepcValidVecComp(w); \
792:   for (i=0;i<vs->n->n;i++) { \
793:     __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i]); \
794:   } \
795:   return(0);\
796: }

800: __FUNC_TEMPLATE3__(Copy)

804: __FUNC_TEMPLATE3__(Swap)


808: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w,Vec z) \
809: { \
810:   PetscErrorCode  ierr; \
811:   Vec_Comp        *vs = (Vec_Comp*)v->data, \
812:                   *ws = (Vec_Comp*)w->data, \
813:                   *zs = (Vec_Comp*)z->data; \
814:   PetscInt        i; \
815: \
817:   SlepcValidVecComp(v); \
818:   SlepcValidVecComp(w); \
819:   SlepcValidVecComp(z); \
820:   for (i=0;i<vs->n->n;i++) { \
821:     __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i],zs->x[i]); \
822:   } \
823:   return(0);\
824: }

828: __FUNC_TEMPLATE4__(PointwiseMax)

832: __FUNC_TEMPLATE4__(PointwiseMaxAbs)

836: __FUNC_TEMPLATE4__(PointwiseMin)

840: __FUNC_TEMPLATE4__(PointwiseMult)

844: __FUNC_TEMPLATE4__(PointwiseDivide)