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)