Actual source code: veccomp.c

slepc-3.17.0 2022-03-31
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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>

 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: static inline void SumNorm2(PetscReal*,PetscReal*,PetscReal*,PetscReal*);
 20: static inline PetscReal GetNorm2(PetscReal,PetscReal);
 21: 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: 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: static inline PetscReal GetNorm2(PetscReal ssq,PetscReal scale)
 44: {
 45:   return scale*PetscSqrtReal(ssq);
 46: }

 48: 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 MPIAPI 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;

 69:   if (*datatype == MPIU_NORM2) {
 70:     for (i=0;i<count;i++) {
 71:       SumNorm2(&xin[i*2],&xin[i*2+1],&xout[i*2],&xout[i*2+1]);
 72:     }
 73:   } else if (*datatype == MPIU_NORM1_AND_2) {
 74:     for (i=0;i<count;i++) {
 75:       xout[i*3] += xin[i*3];
 76:       SumNorm2(&xin[i*3+1],&xin[i*3+2],&xout[i*3+1],&xout[i*3+2]);
 77:     }
 78:   } else {
 79:     (*PetscErrorPrintf)("Can only handle MPIU_NORM* data types");
 80:     MPI_Abort(MPI_COMM_WORLD,1);
 81:   }
 82:   return;
 83: }

 85: static PetscErrorCode VecCompNormEnd(void)
 86: {
 87:   MPI_Type_free(&MPIU_NORM2);
 88:   MPI_Type_free(&MPIU_NORM1_AND_2);
 89:   MPI_Op_free(&MPIU_NORM2_SUM);
 90:   VecCompInitialized = PETSC_FALSE;
 91:   PetscFunctionReturn(0);
 92: }

 94: static PetscErrorCode VecCompNormInit(void)
 95: {
 96:   MPI_Type_contiguous(2,MPIU_REAL,&MPIU_NORM2);
 97:   MPI_Type_commit(&MPIU_NORM2);
 98:   MPI_Type_contiguous(3,MPIU_REAL,&MPIU_NORM1_AND_2);
 99:   MPI_Type_commit(&MPIU_NORM1_AND_2);
100:   MPI_Op_create(SlepcSumNorm2_Local,PETSC_TRUE,&MPIU_NORM2_SUM);
101:   PetscRegisterFinalize(VecCompNormEnd);
102:   PetscFunctionReturn(0);
103: }

105: PetscErrorCode VecDestroy_Comp(Vec v)
106: {
107:   Vec_Comp       *vs = (Vec_Comp*)v->data;
108:   PetscInt       i;

110: #if defined(PETSC_USE_LOG)
111:   PetscLogObjectState((PetscObject)v,"Length=%" PetscInt_FMT,v->map->n);
112: #endif
113:   for (i=0;i<vs->nx;i++) VecDestroy(&vs->x[i]);
114:   if (--vs->n->friends <= 0) PetscFree(vs->n);
115:   PetscFree(vs->x);
116:   PetscFree(vs);
117:   PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",NULL);
118:   PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",NULL);
119:   PetscFunctionReturn(0);
120: }

122: static struct _VecOps DvOps = {
123:   PetscDesignatedInitializer(duplicate,VecDuplicate_Comp),
124:   PetscDesignatedInitializer(duplicatevecs,VecDuplicateVecs_Comp),
125:   PetscDesignatedInitializer(destroyvecs,VecDestroyVecs_Comp),
126:   PetscDesignatedInitializer(dot,VecDot_Comp_MPI),
127:   PetscDesignatedInitializer(mdot,VecMDot_Comp_MPI),
128:   PetscDesignatedInitializer(norm,VecNorm_Comp_MPI),
129:   PetscDesignatedInitializer(tdot,VecTDot_Comp_MPI),
130:   PetscDesignatedInitializer(mtdot,VecMTDot_Comp_MPI),
131:   PetscDesignatedInitializer(scale,VecScale_Comp),
132:   PetscDesignatedInitializer(copy,VecCopy_Comp),
133:   PetscDesignatedInitializer(set,VecSet_Comp),
134:   PetscDesignatedInitializer(swap,VecSwap_Comp),
135:   PetscDesignatedInitializer(axpy,VecAXPY_Comp),
136:   PetscDesignatedInitializer(axpby,VecAXPBY_Comp),
137:   PetscDesignatedInitializer(maxpy,VecMAXPY_Comp),
138:   PetscDesignatedInitializer(aypx,VecAYPX_Comp),
139:   PetscDesignatedInitializer(waxpy,VecWAXPY_Comp),
140:   PetscDesignatedInitializer(axpbypcz,VecAXPBYPCZ_Comp),
141:   PetscDesignatedInitializer(pointwisemult,VecPointwiseMult_Comp),
142:   PetscDesignatedInitializer(pointwisedivide,VecPointwiseDivide_Comp),
143:   PetscDesignatedInitializer(setvalues,NULL),
144:   PetscDesignatedInitializer(assemblybegin,NULL),
145:   PetscDesignatedInitializer(assemblyend,NULL),
146:   PetscDesignatedInitializer(getarray,NULL),
147:   PetscDesignatedInitializer(getsize,VecGetSize_Comp),
148:   PetscDesignatedInitializer(getlocalsize,VecGetLocalSize_Comp),
149:   PetscDesignatedInitializer(restorearray,NULL),
150:   PetscDesignatedInitializer(max,VecMax_Comp),
151:   PetscDesignatedInitializer(min,VecMin_Comp),
152:   PetscDesignatedInitializer(setrandom,VecSetRandom_Comp),
153:   PetscDesignatedInitializer(setoption,NULL),
154:   PetscDesignatedInitializer(setvaluesblocked,NULL),
155:   PetscDesignatedInitializer(destroy,VecDestroy_Comp),
156:   PetscDesignatedInitializer(view,VecView_Comp),
157:   PetscDesignatedInitializer(placearray,NULL),
158:   PetscDesignatedInitializer(replacearray,NULL),
159:   PetscDesignatedInitializer(dot_local,VecDot_Comp_Seq),
160:   PetscDesignatedInitializer(tdot_local,VecTDot_Comp_Seq),
161:   PetscDesignatedInitializer(norm_local,VecNorm_Comp_Seq),
162:   PetscDesignatedInitializer(mdot_local,VecMDot_Comp_Seq),
163:   PetscDesignatedInitializer(mtdot_local,VecMTDot_Comp_Seq),
164:   PetscDesignatedInitializer(load,NULL),
165:   PetscDesignatedInitializer(reciprocal,VecReciprocal_Comp),
166:   PetscDesignatedInitializer(conjugate,VecConjugate_Comp),
167:   PetscDesignatedInitializer(setlocaltoglobalmapping,NULL),
168:   PetscDesignatedInitializer(setvalueslocal,NULL),
169:   PetscDesignatedInitializer(resetarray,NULL),
170:   PetscDesignatedInitializer(setfromoptions,NULL),
171:   PetscDesignatedInitializer(maxpointwisedivide,VecMaxPointwiseDivide_Comp),
172:   PetscDesignatedInitializer(pointwisemax,VecPointwiseMax_Comp),
173:   PetscDesignatedInitializer(pointwisemaxabs,VecPointwiseMaxAbs_Comp),
174:   PetscDesignatedInitializer(pointwisemin,VecPointwiseMin_Comp),
175:   PetscDesignatedInitializer(getvalues,NULL),
176:   PetscDesignatedInitializer(sqrt,VecSqrtAbs_Comp),
177:   PetscDesignatedInitializer(abs,VecAbs_Comp),
178:   PetscDesignatedInitializer(exp,VecExp_Comp),
179:   PetscDesignatedInitializer(log,VecLog_Comp),
180:   PetscDesignatedInitializer(shift,VecShift_Comp),
181:   PetscDesignatedInitializer(create,NULL),
182:   PetscDesignatedInitializer(stridegather,NULL),
183:   PetscDesignatedInitializer(stridescatter,NULL),
184:   PetscDesignatedInitializer(dotnorm2,VecDotNorm2_Comp_MPI),
185:   PetscDesignatedInitializer(getsubvector,NULL),
186:   PetscDesignatedInitializer(restoresubvector,NULL),
187:   PetscDesignatedInitializer(getarrayread,NULL),
188:   PetscDesignatedInitializer(restorearrayread,NULL),
189:   PetscDesignatedInitializer(stridesubsetgather,NULL),
190:   PetscDesignatedInitializer(stridesubsetscatter,NULL),
191:   PetscDesignatedInitializer(viewnative,NULL),
192:   PetscDesignatedInitializer(loadnative,NULL),
193:   PetscDesignatedInitializer(getlocalvector,NULL)
194: };

196: PetscErrorCode VecDuplicateVecs_Comp(Vec w,PetscInt m,Vec *V[])
197: {
198:   PetscInt       i;

203:   PetscMalloc1(m,V);
204:   for (i=0;i<m;i++) VecDuplicate(w,*V+i);
205:   PetscFunctionReturn(0);
206: }

208: PetscErrorCode VecDestroyVecs_Comp(PetscInt m,Vec v[])
209: {
210:   PetscInt       i;

214:   for (i=0;i<m;i++) VecDestroy(&v[i]);
215:   PetscFree(v);
216:   PetscFunctionReturn(0);
217: }

219: static PetscErrorCode VecCreate_Comp_Private(Vec v,Vec *x,PetscInt nx,PetscBool x_to_me,Vec_Comp_N *n)
220: {
221:   Vec_Comp       *s;
222:   PetscInt       N=0,lN=0,i,k;

224:   if (!VecCompInitialized) {
225:     VecCompInitialized = PETSC_TRUE;
226:     VecRegister(VECCOMP,VecCreate_Comp);
227:     VecCompNormInit();
228:   }

230:   /* Allocate a new Vec_Comp */
231:   if (v->data) PetscFree(v->data);
232:   PetscNewLog(v,&s);
233:   PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
234:   v->data        = (void*)s;
235:   v->petscnative = PETSC_FALSE;

237:   /* Allocate the array of Vec, if it is needed to be done */
238:   if (!x_to_me) {
239:     if (nx) PetscMalloc1(nx,&s->x);
240:     if (x) PetscArraycpy(s->x,x,nx);
241:   } else s->x = x;

243:   s->nx = nx;

245:   if (nx && x) {
246:     /* Allocate the shared structure, if it is not given */
247:     if (!n) {
248:       for (i=0;i<nx;i++) {
249:         VecGetSize(x[i],&k);
250:         N+= k;
251:         VecGetLocalSize(x[i],&k);
252:         lN+= k;
253:       }
254:       PetscNewLog(v,&n);
255:       s->n = n;
256:       n->n = nx;
257:       n->N = N;
258:       n->lN = lN;
259:       n->friends = 1;
260:     } else { /* If not, check in the vector in the shared structure */
261:       s->n = n;
262:       s->n->friends++;
263:     }

265:     /* Set the virtual sizes as the real sizes of the vector */
266:     VecSetSizes(v,s->n->lN,s->n->N);
267:   }

269:   PetscObjectChangeTypeName((PetscObject)v,VECCOMP);
270:   PetscObjectComposeFunction((PetscObject)v,"VecCompSetSubVecs_C",VecCompSetSubVecs_Comp);
271:   PetscObjectComposeFunction((PetscObject)v,"VecCompGetSubVecs_C",VecCompGetSubVecs_Comp);
272:   PetscFunctionReturn(0);
273: }

275: SLEPC_EXTERN PetscErrorCode VecCreate_Comp(Vec V)
276: {
277:   VecCreate_Comp_Private(V,NULL,0,PETSC_FALSE,NULL);
278:   PetscFunctionReturn(0);
279: }

281: /*@C
282:    VecCreateComp - Creates a new vector containing several subvectors,
283:    each stored separately.

285:    Collective

287:    Input Parameters:
288: +  comm - communicator for the new Vec
289: .  Nx   - array of (initial) global sizes of child vectors
290: .  n    - number of child vectors
291: .  t    - type of the child vectors
292: -  Vparent - (optional) template vector

294:    Output Parameter:
295: .  V - new vector

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

301:    Level: developer

303: .seealso: VecCreateCompWithVecs(), VecCompSetSubVecs()
304: @*/
305: PetscErrorCode VecCreateComp(MPI_Comm comm,PetscInt *Nx,PetscInt n,VecType t,Vec Vparent,Vec *V)
306: {
307:   Vec            *x;
308:   PetscInt       i;

310:   VecCreate(comm,V);
311:   PetscMalloc1(n,&x);
312:   PetscLogObjectMemory((PetscObject)*V,n*sizeof(Vec));
313:   for (i=0;i<n;i++) {
314:     VecCreate(comm,&x[i]);
315:     VecSetSizes(x[i],PETSC_DECIDE,Nx[i]);
316:     VecSetType(x[i],t);
317:   }
318:   VecCreate_Comp_Private(*V,x,n,PETSC_TRUE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
319:   PetscFunctionReturn(0);
320: }

322: /*@C
323:    VecCreateCompWithVecs - Creates a new vector containing several subvectors,
324:    each stored separately, from an array of Vecs.

326:    Collective on x

328:    Input Parameters:
329: +  x - array of Vecs
330: .  n - number of child vectors
331: -  Vparent - (optional) template vector

333:    Output Parameter:
334: .  V - new vector

336:    Level: developer

338: .seealso: VecCreateComp()
339: @*/
340: PetscErrorCode VecCreateCompWithVecs(Vec *x,PetscInt n,Vec Vparent,Vec *V)
341: {
342:   PetscInt       i;

347:   VecCreate(PetscObjectComm((PetscObject)x[0]),V);
348:   for (i=0;i<n;i++) PetscObjectReference((PetscObject)x[i]);
349:   VecCreate_Comp_Private(*V,x,n,PETSC_FALSE,Vparent?((Vec_Comp*)Vparent->data)->n:NULL);
350:   PetscFunctionReturn(0);
351: }

353: PetscErrorCode VecDuplicate_Comp(Vec win,Vec *V)
354: {
355:   Vec            *x;
356:   PetscInt       i;
357:   Vec_Comp       *s = (Vec_Comp*)win->data;

359:   SlepcValidVecComp(win,1);
360:   VecCreate(PetscObjectComm((PetscObject)win),V);
361:   PetscMalloc1(s->nx,&x);
362:   PetscLogObjectMemory((PetscObject)*V,s->nx*sizeof(Vec));
363:   for (i=0;i<s->nx;i++) {
364:     if (s->x[i]) VecDuplicate(s->x[i],&x[i]);
365:     else x[i] = NULL;
366:   }
367:   VecCreate_Comp_Private(*V,x,s->nx,PETSC_TRUE,s->n);
368:   PetscFunctionReturn(0);
369: }

371: static PetscErrorCode VecCompGetSubVecs_Comp(Vec win,PetscInt *n,const Vec **x)
372: {
373:   Vec_Comp *s = (Vec_Comp*)win->data;

375:   if (x) *x = s->x;
376:   if (n) *n = s->n->n;
377:   PetscFunctionReturn(0);
378: }

380: /*@C
381:    VecCompGetSubVecs - Returns the entire array of vectors defining a
382:    compound vector.

384:    Collective on win

386:    Input Parameter:
387: .  win - compound vector

389:    Output Parameters:
390: +  n - number of child vectors
391: -  x - array of child vectors

393:    Level: developer

395: .seealso: VecCreateComp()
396: @*/
397: PetscErrorCode VecCompGetSubVecs(Vec win,PetscInt *n,const Vec **x)
398: {
400:   PetscUseMethod(win,"VecCompGetSubVecs_C",(Vec,PetscInt*,const Vec**),(win,n,x));
401:   PetscFunctionReturn(0);
402: }

404: static PetscErrorCode VecCompSetSubVecs_Comp(Vec win,PetscInt n,Vec *x)
405: {
406:   Vec_Comp       *s = (Vec_Comp*)win->data;
407:   PetscInt       i,N,nlocal;
408:   Vec_Comp_N     *nn;

411:   if (!s->nx) {
412:     /* vector has been created via VecCreate+VecSetType+VecSetSizes, so allocate data structures */
413:     PetscMalloc1(n,&s->x);
414:     PetscLogObjectMemory((PetscObject)win,n*sizeof(Vec));
415:     VecGetSize(win,&N);
417:     VecGetLocalSize(win,&nlocal);
419:     s->nx = n;
420:     for (i=0;i<n;i++) {
421:       VecCreate(PetscObjectComm((PetscObject)win),&s->x[i]);
422:       VecSetSizes(s->x[i],nlocal/n,N/n);
423:       VecSetFromOptions(s->x[i]);
424:     }
425:     if (!s->n) {
426:       PetscNewLog(win,&nn);
427:       s->n = nn;
428:       nn->N = N;
429:       nn->lN = nlocal;
430:       nn->friends = 1;
431:     }
433:   if (x) PetscArraycpy(s->x,x,n);
434:   s->n->n = n;
435:   PetscFunctionReturn(0);
436: }

438: /*@C
439:    VecCompSetSubVecs - Resets the number of subvectors defining a compound vector,
440:    or replaces the subvectors.

442:    Collective on win

444:    Input Parameters:
445: +  win - compound vector
446: .  n - number of child vectors
447: -  x - array of child vectors

449:    Note:
450:    It is not possible to increase the number of subvectors with respect to the
451:    number set at its creation.

453:    Level: developer

455: .seealso: VecCreateComp(), VecCompGetSubVecs()
456: @*/
457: PetscErrorCode VecCompSetSubVecs(Vec win,PetscInt n,Vec *x)
458: {
461:   PetscTryMethod(win,"VecCompSetSubVecs_C",(Vec,PetscInt,Vec*),(win,n,x));
462:   PetscFunctionReturn(0);
463: }

465: PetscErrorCode VecAXPY_Comp(Vec v,PetscScalar alpha,Vec w)
466: {
467:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
468:   PetscInt       i;

470:   SlepcValidVecComp(v,1);
471:   SlepcValidVecComp(w,3);
472:   for (i=0;i<vs->n->n;i++) VecAXPY(vs->x[i],alpha,ws->x[i]);
473:   PetscFunctionReturn(0);
474: }

476: PetscErrorCode VecAYPX_Comp(Vec v,PetscScalar alpha,Vec w)
477: {
478:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
479:   PetscInt       i;

481:   SlepcValidVecComp(v,1);
482:   SlepcValidVecComp(w,3);
483:   for (i=0;i<vs->n->n;i++) VecAYPX(vs->x[i],alpha,ws->x[i]);
484:   PetscFunctionReturn(0);
485: }

487: PetscErrorCode VecAXPBY_Comp(Vec v,PetscScalar alpha,PetscScalar beta,Vec w)
488: {
489:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
490:   PetscInt       i;

492:   SlepcValidVecComp(v,1);
493:   SlepcValidVecComp(w,4);
494:   for (i=0;i<vs->n->n;i++) VecAXPBY(vs->x[i],alpha,beta,ws->x[i]);
495:   PetscFunctionReturn(0);
496: }

498: PetscErrorCode VecMAXPY_Comp(Vec v,PetscInt n,const PetscScalar *alpha,Vec *w)
499: {
500:   Vec_Comp       *vs = (Vec_Comp*)v->data;
501:   Vec            *wx;
502:   PetscInt       i,j;

504:   SlepcValidVecComp(v,1);
505:   for (i=0;i<n;i++) SlepcValidVecComp(w[i],4);

507:   PetscMalloc1(n,&wx);

509:   for (j=0;j<vs->n->n;j++) {
510:     for (i=0;i<n;i++) wx[i] = ((Vec_Comp*)w[i]->data)->x[j];
511:     VecMAXPY(vs->x[j],n,alpha,wx);
512:   }

514:   PetscFree(wx);
515:   PetscFunctionReturn(0);
516: }

518: PetscErrorCode VecWAXPY_Comp(Vec v,PetscScalar alpha,Vec w,Vec z)
519: {
520:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
521:   PetscInt       i;

523:   SlepcValidVecComp(v,1);
524:   SlepcValidVecComp(w,3);
525:   SlepcValidVecComp(z,4);
526:   for (i=0;i<vs->n->n;i++) VecWAXPY(vs->x[i],alpha,ws->x[i],zs->x[i]);
527:   PetscFunctionReturn(0);
528: }

530: PetscErrorCode VecAXPBYPCZ_Comp(Vec v,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec w,Vec z)
531: {
532:   Vec_Comp        *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data,*zs = (Vec_Comp*)z->data;
533:   PetscInt        i;

535:   SlepcValidVecComp(v,1);
536:   SlepcValidVecComp(w,5);
537:   SlepcValidVecComp(z,6);
538:   for (i=0;i<vs->n->n;i++) VecAXPBYPCZ(vs->x[i],alpha,beta,gamma,ws->x[i],zs->x[i]);
539:   PetscFunctionReturn(0);
540: }

542: PetscErrorCode VecGetSize_Comp(Vec v,PetscInt *size)
543: {
544:   Vec_Comp *vs = (Vec_Comp*)v->data;

547:   if (vs->n) {
548:     SlepcValidVecComp(v,1);
549:     *size = vs->n->N;
550:   } else *size = v->map->N;
551:   PetscFunctionReturn(0);
552: }

554: PetscErrorCode VecGetLocalSize_Comp(Vec v,PetscInt *size)
555: {
556:   Vec_Comp *vs = (Vec_Comp*)v->data;

559:   if (vs->n) {
560:     SlepcValidVecComp(v,1);
561:     *size = vs->n->lN;
562:   } else *size = v->map->n;
563:   PetscFunctionReturn(0);
564: }

566: PetscErrorCode VecMax_Comp(Vec v,PetscInt *idx,PetscReal *z)
567: {
568:   Vec_Comp       *vs = (Vec_Comp*)v->data;
569:   PetscInt       idxp,s=0,s0;
570:   PetscReal      zp,z0;
571:   PetscInt       i;

573:   SlepcValidVecComp(v,1);
574:   if (!idx && !z) PetscFunctionReturn(0);

576:   if (vs->n->n > 0) VecMax(vs->x[0],idx?&idxp:NULL,&zp);
577:   else {
578:     zp = PETSC_MIN_REAL;
579:     if (idx) idxp = -1;
580:   }
581:   for (i=1;i<vs->n->n;i++) {
582:     VecGetSize(vs->x[i-1],&s0);
583:     s += s0;
584:     VecMax(vs->x[i],idx?&idxp:NULL,&z0);
585:     if (zp < z0) {
586:       if (idx) *idx = s+idxp;
587:       zp = z0;
588:     }
589:   }
590:   if (z) *z = zp;
591:   PetscFunctionReturn(0);
592: }

594: PetscErrorCode VecMin_Comp(Vec v,PetscInt *idx,PetscReal *z)
595: {
596:   Vec_Comp       *vs = (Vec_Comp*)v->data;
597:   PetscInt       idxp,s=0,s0;
598:   PetscReal      zp,z0;
599:   PetscInt       i;

601:   SlepcValidVecComp(v,1);
602:   if (!idx && !z) PetscFunctionReturn(0);

604:   if (vs->n->n > 0) VecMin(vs->x[0],idx?&idxp:NULL,&zp);
605:   else {
606:     zp = PETSC_MAX_REAL;
607:     if (idx) idxp = -1;
608:   }
609:   for (i=1;i<vs->n->n;i++) {
610:     VecGetSize(vs->x[i-1],&s0);
611:     s += s0;
612:     VecMin(vs->x[i],idx?&idxp:NULL,&z0);
613:     if (zp > z0) {
614:       if (idx) *idx = s+idxp;
615:       zp = z0;
616:     }
617:   }
618:   if (z) *z = zp;
619:   PetscFunctionReturn(0);
620: }

622: PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v,Vec w,PetscReal *m)
623: {
624:   Vec_Comp       *vs = (Vec_Comp*)v->data,*ws = (Vec_Comp*)w->data;
625:   PetscReal      work;
626:   PetscInt       i;

628:   SlepcValidVecComp(v,1);
629:   SlepcValidVecComp(w,2);
630:   if (!m || vs->n->n == 0) PetscFunctionReturn(0);
631:   VecMaxPointwiseDivide(vs->x[0],ws->x[0],m);
632:   for (i=1;i<vs->n->n;i++) {
633:     VecMaxPointwiseDivide(vs->x[i],ws->x[i],&work);
634:     *m = PetscMax(*m,work);
635:   }
636:   PetscFunctionReturn(0);
637: }


644: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v) \
645: { \
646:   Vec_Comp        *vs = (Vec_Comp*)v->data; \
647:   PetscInt        i; \
648: \
649:   SlepcValidVecComp(v,1); \
650:   for (i=0;i<vs->n->n;i++) { \
651:     __COMPOSE2__(Vec,NAME)(vs->x[i]); \
652:   } \
653:   PetscFunctionReturn(0);\
654: }

656: __FUNC_TEMPLATE1__(Conjugate)
657: __FUNC_TEMPLATE1__(Reciprocal)
658: __FUNC_TEMPLATE1__(SqrtAbs)
659: __FUNC_TEMPLATE1__(Abs)
660: __FUNC_TEMPLATE1__(Exp)
661: __FUNC_TEMPLATE1__(Log)

664: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,T0 __a) \
665: { \
666:   Vec_Comp        *vs = (Vec_Comp*)v->data; \
667:   PetscInt        i; \
668: \
669:   SlepcValidVecComp(v,1); \
670:   for (i=0;i<vs->n->n;i++) { \
671:     __COMPOSE2__(Vec,NAME)(vs->x[i],__a); \
672:   } \
673:   PetscFunctionReturn(0);\
674: }

676: __FUNC_TEMPLATE2__(Set,PetscScalar)
677: __FUNC_TEMPLATE2__(View,PetscViewer)
678: __FUNC_TEMPLATE2__(Scale,PetscScalar)
679: __FUNC_TEMPLATE2__(SetRandom,PetscRandom)
680: __FUNC_TEMPLATE2__(Shift,PetscScalar)

683: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w) \
684: { \
685:   Vec_Comp        *vs = (Vec_Comp*)v->data,\
686:                   *ws = (Vec_Comp*)w->data; \
687:   PetscInt        i; \
688: \
689:   SlepcValidVecComp(v,1); \
690:   SlepcValidVecComp(w,2); \
691:   for (i=0;i<vs->n->n;i++) { \
692:     __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i]); \
693:   } \
694:   PetscFunctionReturn(0);\
695: }

697: __FUNC_TEMPLATE3__(Copy)
698: __FUNC_TEMPLATE3__(Swap)

701: PetscErrorCode __COMPOSE3__(Vec,NAME,_Comp)(Vec v,Vec w,Vec z) \
702: { \
703:   Vec_Comp        *vs = (Vec_Comp*)v->data, \
704:                   *ws = (Vec_Comp*)w->data, \
705:                   *zs = (Vec_Comp*)z->data; \
706:   PetscInt        i; \
707: \
708:   SlepcValidVecComp(v,1); \
709:   SlepcValidVecComp(w,2); \
710:   SlepcValidVecComp(z,3); \
711:   for (i=0;i<vs->n->n;i++) { \
712:     __COMPOSE2__(Vec,NAME)(vs->x[i],ws->x[i],zs->x[i]); \
713:   } \
714:   PetscFunctionReturn(0);\
715: }

717: __FUNC_TEMPLATE4__(PointwiseMax)
718: __FUNC_TEMPLATE4__(PointwiseMaxAbs)
719: __FUNC_TEMPLATE4__(PointwiseMin)
720: __FUNC_TEMPLATE4__(PointwiseMult)
721: __FUNC_TEMPLATE4__(PointwiseDivide)