Actual source code: bvec2.c

  1: #define PETSCVEC_DLL
  2: /*
  3:    Implements the sequential vectors.
  4: */

 6:  #include vecimpl.h
 7:  #include src/vec/impls/dvecimpl.h
 8:  #include src/inline/dot.h
 9:  #include petscblaslapack.h
 10: #if defined(PETSC_HAVE_PNETCDF)
 12: #include "pnetcdf.h"
 14: #endif

 18: PetscErrorCode VecNorm_Seq(Vec xin,NormType type,PetscReal* z)
 19: {
 20:   PetscScalar    *xx;
 22:   PetscInt       n = xin->n;
 23:   PetscBLASInt   bn = (PetscBLASInt)n,one = 1;

 26:   if (type == NORM_2 || type == NORM_FROBENIUS) {
 27:     VecGetArray(xin,&xx);
 28:     /*
 29:       This is because the Fortran BLAS 1 Norm is very slow! 
 30:     */
 31: #if defined(PETSC_HAVE_SLOW_BLAS_NORM2)
 32: #if defined(PETSC_USE_FORTRAN_KERNEL_NORM)
 33:     fortrannormsqr_(xx,&n,z);
 34:     *z = sqrt(*z);
 35: #elif defined(PETSC_USE_UNROLLED_NORM)
 36:     {
 37:     PetscReal work = 0.0;
 38:     switch (n & 0x3) {
 39:       case 3: work += PetscRealPart(xx[0]*PetscConj(xx[0])); xx++;
 40:       case 2: work += PetscRealPart(xx[0]*PetscConj(xx[0])); xx++;
 41:       case 1: work += PetscRealPart(xx[0]*PetscConj(xx[0])); xx++; n -= 4;
 42:     }
 43:     while (n>0) {
 44:       work += PetscRealPart(xx[0]*PetscConj(xx[0])+xx[1]*PetscConj(xx[1])+
 45:                         xx[2]*PetscConj(xx[2])+xx[3]*PetscConj(xx[3]));
 46:       xx += 4; n -= 4;
 47:     }
 48:     *z = sqrt(work);}
 49: #else
 50:     {
 51:       PetscInt         i;
 52:       PetscScalar sum=0.0;
 53:       for (i=0; i<n; i++) {
 54:         sum += (xx[i])*(PetscConj(xx[i]));
 55:       }
 56:       *z = sqrt(PetscRealPart(sum));
 57:     }
 58: #endif
 59: #else
 60:     *z = BLASnrm2_(&bn,xx,&one);
 61: #endif
 62:     VecRestoreArray(xin,&xx);
 63:     PetscLogFlops(2*n-1);
 64:   } else if (type == NORM_INFINITY) {
 65:     PetscInt          i;
 66:     PetscReal    max = 0.0,tmp;

 68:     VecGetArray(xin,&xx);
 69:     for (i=0; i<n; i++) {
 70:       if ((tmp = PetscAbsScalar(*xx)) > max) max = tmp;
 71:       /* check special case of tmp == NaN */
 72:       if (tmp != tmp) {max = tmp; break;}
 73:       xx++;
 74:     }
 75:     VecRestoreArray(xin,&xx);
 76:     *z   = max;
 77:   } else if (type == NORM_1) {
 78:     VecGetArray(xin,&xx);
 79:     *z = BLASasum_(&bn,xx,&one);
 80:     VecRestoreArray(xin,&xx);
 81:     PetscLogFlops(n-1);
 82:   } else if (type == NORM_1_AND_2) {
 83:     VecNorm_Seq(xin,NORM_1,z);
 84:     VecNorm_Seq(xin,NORM_2,z+1);
 85:   }
 86:   return(0);
 87: }

 91: PetscErrorCode VecView_Seq_File(Vec xin,PetscViewer viewer)
 92: {
 93:   Vec_Seq           *x = (Vec_Seq *)xin->data;
 94:   PetscErrorCode    ierr;
 95:   PetscInt          i,n = xin->n;
 96:   const char        *name;
 97:   PetscViewerFormat format;

100:   PetscViewerGetFormat(viewer,&format);
101:   if (format == PETSC_VIEWER_ASCII_MATLAB) {
102:     PetscObjectGetName((PetscObject)xin,&name);
103:     PetscViewerASCIIPrintf(viewer,"%s = [\n",name);
104:     for (i=0; i<n; i++) {
105: #if defined(PETSC_USE_COMPLEX)
106:       if (PetscImaginaryPart(x->array[i]) > 0.0) {
107:         PetscViewerASCIIPrintf(viewer,"%18.16e + %18.16ei\n",PetscRealPart(x->array[i]),PetscImaginaryPart(x->array[i]));
108:       } else if (PetscImaginaryPart(x->array[i]) < 0.0) {
109:         PetscViewerASCIIPrintf(viewer,"%18.16e - %18.16ei\n",PetscRealPart(x->array[i]),-PetscImaginaryPart(x->array[i]));
110:       } else {
111:         PetscViewerASCIIPrintf(viewer,"%18.16e\n",PetscRealPart(x->array[i]));
112:       }
113: #else
114:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",x->array[i]);
115: #endif
116:     }
117:     PetscViewerASCIIPrintf(viewer,"];\n");
118:   } else if (format == PETSC_VIEWER_ASCII_SYMMODU) {
119:     for (i=0; i<n; i++) {
120: #if defined(PETSC_USE_COMPLEX)
121:       PetscViewerASCIIPrintf(viewer,"%18.16e %18.16e\n",PetscRealPart(x->array[i]),PetscImaginaryPart(x->array[i]));
122: #else
123:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",x->array[i]);
124: #endif
125:     }
126:   } else {
127:     for (i=0; i<n; i++) {
128:       if (format == PETSC_VIEWER_ASCII_INDEX) {
129:         PetscViewerASCIIPrintf(viewer,"%D: ",i);
130:       }
131: #if defined(PETSC_USE_COMPLEX)
132:       if (PetscImaginaryPart(x->array[i]) > 0.0) {
133:         PetscViewerASCIIPrintf(viewer,"%g + %g i\n",PetscRealPart(x->array[i]),PetscImaginaryPart(x->array[i]));
134:       } else if (PetscImaginaryPart(x->array[i]) < 0.0) {
135:         PetscViewerASCIIPrintf(viewer,"%g - %g i\n",PetscRealPart(x->array[i]),-PetscImaginaryPart(x->array[i]));
136:       } else {
137:         PetscViewerASCIIPrintf(viewer,"%g\n",PetscRealPart(x->array[i]));
138:       }
139: #else
140:       PetscViewerASCIIPrintf(viewer,"%g\n",x->array[i]);
141: #endif
142:     }
143:   }
144:   PetscViewerFlush(viewer);
145:   return(0);
146: }

150: static PetscErrorCode VecView_Seq_Draw_LG(Vec xin,PetscViewer v)
151: {
152:   Vec_Seq        *x = (Vec_Seq *)xin->data;
154:   PetscInt       i,n = xin->n;
155:   PetscDraw      win;
156:   PetscReal      *xx;
157:   PetscDrawLG    lg;

160:   PetscViewerDrawGetDrawLG(v,0,&lg);
161:   PetscDrawLGGetDraw(lg,&win);
162:   PetscDrawCheckResizedWindow(win);
163:   PetscDrawLGReset(lg);
164:   PetscMalloc((n+1)*sizeof(PetscReal),&xx);
165:   for (i=0; i<n; i++) {
166:     xx[i] = (PetscReal) i;
167:   }
168: #if !defined(PETSC_USE_COMPLEX)
169:   PetscDrawLGAddPoints(lg,n,&xx,&x->array);
170: #else 
171:   {
172:     PetscReal *yy;
173:     PetscMalloc((n+1)*sizeof(PetscReal),&yy);
174:     for (i=0; i<n; i++) {
175:       yy[i] = PetscRealPart(x->array[i]);
176:     }
177:     PetscDrawLGAddPoints(lg,n,&xx,&yy);
178:     PetscFree(yy);
179:   }
180: #endif
181:   PetscFree(xx);
182:   PetscDrawLGDraw(lg);
183:   PetscDrawSynchronizedFlush(win);
184:   return(0);
185: }

189: static PetscErrorCode VecView_Seq_Draw(Vec xin,PetscViewer v)
190: {
191:   PetscErrorCode    ierr;
192:   PetscDraw         draw;
193:   PetscTruth        isnull;
194:   PetscViewerFormat format;

197:   PetscViewerDrawGetDraw(v,0,&draw);
198:   PetscDrawIsNull(draw,&isnull); if (isnull) return(0);
199: 
200:   PetscViewerGetFormat(v,&format);
201:   /*
202:      Currently it only supports drawing to a line graph */
203:   if (format != PETSC_VIEWER_DRAW_LG) {
204:     PetscViewerPushFormat(v,PETSC_VIEWER_DRAW_LG);
205:   }
206:   VecView_Seq_Draw_LG(xin,v);
207:   if (format != PETSC_VIEWER_DRAW_LG) {
208:     PetscViewerPopFormat(v);
209:   }

211:   return(0);
212: }

216: static PetscErrorCode VecView_Seq_Binary(Vec xin,PetscViewer viewer)
217: {
218:   Vec_Seq        *x = (Vec_Seq *)xin->data;
220:   int            fdes;
221:   PetscInt       n = xin->n,cookie=VEC_FILE_COOKIE;
222:   FILE           *file;

225:   PetscViewerBinaryGetDescriptor(viewer,&fdes);
226:   /* Write vector header */
227:   PetscBinaryWrite(fdes,&cookie,1,PETSC_INT,PETSC_FALSE);
228:   PetscBinaryWrite(fdes,&n,1,PETSC_INT,PETSC_FALSE);

230:   /* Write vector contents */
231:   PetscBinaryWrite(fdes,x->array,n,PETSC_SCALAR,PETSC_FALSE);

233:   PetscViewerBinaryGetInfoPointer(viewer,&file);
234:   if (file && xin->bs > 1) {
235:     if (xin->prefix) {
236:       PetscFPrintf(PETSC_COMM_SELF,file,"-%s_vecload_block_size %D\n",xin->prefix,xin->bs);
237:     } else {
238:       PetscFPrintf(PETSC_COMM_SELF,file,"-vecload_block_size %D\n",xin->bs);
239:     }
240:   }
241:   return(0);
242: }

244: #if defined(PETSC_HAVE_PNETCDF)
247: PetscErrorCode VecView_Seq_Netcdf(Vec xin,PetscViewer v)
248: {
250:   int            n = xin->n,ncid,xdim,xdim_num=1,xin_id,xstart=0;
251:   MPI_Comm       comm = xin->comm;
252:   PetscScalar    *values,*xarray;

255: #if !defined(PETSC_USE_COMPLEX)
256:   VecGetArray(xin,&xarray);
257:   PetscViewerNetcdfGetID(v,&ncid);
258:   if (ncid < 0) SETERRQ(PETSC_ERR_ORDER,"First call PetscViewerNetcdfOpen to create NetCDF dataset");
259:   /* define dimensions */
260:   ncmpi_def_dim(ncid,"PETSc_Vector_Global_Size",n,&xdim);
261:   /* define variables */
262:   ncmpi_def_var(ncid,"PETSc_Vector_Seq",NC_DOUBLE,xdim_num,&xdim,&xin_id);
263:   /* leave define mode */
264:   ncmpi_enddef(ncid);
265:   /* store the vector */
266:   VecGetOwnershipRange(xin,&xstart,PETSC_NULL);
267:   ncmpi_put_vara_double_all(ncid,xin_id,(const size_t*)&xstart,(const size_t*)&n,xarray);
268: #else 
269:     PetscPrintf(PETSC_COMM_WORLD,"NetCDF viewer not supported for complex numbers\n");
270: #endif
271:   return(0);
272: }
273: #endif

275: #if defined(PETSC_HAVE_MATLAB)
276: #include "mat.h"   /* Matlab include file */
280: PetscErrorCode VecView_Seq_Matlab(Vec vec,PetscViewer viewer)
281: {
283:   PetscInt       n;
284:   PetscScalar    *array;
285: 
287:   VecGetLocalSize(vec,&n);
288:   PetscObjectName((PetscObject)vec);
289:   VecGetArray(vec,&array);
290:   PetscViewerMatlabPutArray(viewer,n,1,array,vec->name);
291:   VecRestoreArray(vec,&array);
292:   return(0);
293: }
295: #endif

299: PetscErrorCode VecView_Seq(Vec xin,PetscViewer viewer)
300: {
302:   PetscTruth     isdraw,iascii,issocket,isbinary;
303: #if defined(PETSC_HAVE_MATHEMATICA)
304:   PetscTruth     ismathematica;
305: #endif
306: #if defined(PETSC_HAVE_PNETCDF)
307:   PetscTruth     isnetcdf;
308: #endif
309: #if defined(PETSC_HAVE_MATLAB)
310:   PetscTruth     ismatlab;
311: #endif

314:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_DRAW,&isdraw);
315:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
316:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_SOCKET,&issocket);
317:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_BINARY,&isbinary);
318: #if defined(PETSC_HAVE_MATHEMATICA)
319:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_MATHEMATICA,&ismathematica);
320: #endif
321: #if defined(PETSC_HAVE_PNETCDF)
322:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_NETCDF,&isnetcdf);
323: #endif
324: #if defined(PETSC_HAVE_MATLAB)
325:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_MATLAB,&ismatlab);
326: #endif

328:   if (isdraw){
329:     VecView_Seq_Draw(xin,viewer);
330:   } else if (iascii){
331:     VecView_Seq_File(xin,viewer);
332: #if defined(PETSC_USE_SOCKET_VIEWER)
333:   } else if (issocket) {
334:     Vec_Seq        *x = (Vec_Seq *)xin->data;
335:     PetscViewerSocketPutScalar(viewer,xin->n,1,x->array);
336: #endif
337:   } else if (isbinary) {
338:     VecView_Seq_Binary(xin,viewer);
339: #if defined(PETSC_HAVE_MATHEMATICA)
340:   } else if (ismathematica) {
341:     PetscViewerMathematicaPutVector(viewer,xin);
342: #endif
343: #if defined(PETSC_HAVE_PNETCDF)
344:   } else if (isnetcdf) {
345:     VecView_Seq_Netcdf(xin,viewer);
346: #endif
347: #if defined(PETSC_HAVE_MATLAB)
348:   } else if (ismatlab) {
349:     VecView_Seq_Matlab(xin,viewer);
350: #endif
351:   } else {
352:     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported by this vector object",((PetscObject)viewer)->type_name);
353:   }
354:   return(0);
355: }

359: PetscErrorCode VecGetValues_Seq(Vec xin,PetscInt ni,const PetscInt ix[],PetscScalar y[])
360: {
361:   Vec_Seq     *x = (Vec_Seq *)xin->data;
362:   PetscScalar *xx = x->array;
363:   PetscInt    i;

366:   for (i=0; i<ni; i++) {
367: #if defined(PETSC_USE_DEBUG)
368:     if (ix[i] < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D cannot be negative",ix[i]);
369:     if (ix[i] >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D to large maximum allowed %D",ix[i],xin->n);
370: #endif
371:     y[i] = xx[ix[i]];
372:   }
373:   return(0);
374: }

378: PetscErrorCode VecSetValues_Seq(Vec xin,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode m)
379: {
380:   Vec_Seq     *x = (Vec_Seq *)xin->data;
381:   PetscScalar *xx = x->array;
382:   PetscInt    i;

385:   if (m == INSERT_VALUES) {
386:     for (i=0; i<ni; i++) {
387:       if (ix[i] < 0) continue;
388: #if defined(PETSC_USE_DEBUG)
389:       if (ix[i] >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",ix[i],xin->n);
390: #endif
391:       xx[ix[i]] = y[i];
392:     }
393:   } else {
394:     for (i=0; i<ni; i++) {
395:       if (ix[i] < 0) continue;
396: #if defined(PETSC_USE_DEBUG)
397:       if (ix[i] >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",ix[i],xin->n);
398: #endif
399:       xx[ix[i]] += y[i];
400:     }
401:   }
402:   return(0);
403: }

407: PetscErrorCode VecSetValuesBlocked_Seq(Vec xin,PetscInt ni,const PetscInt ix[],const PetscScalar yin[],InsertMode m)
408: {
409:   Vec_Seq     *x = (Vec_Seq *)xin->data;
410:   PetscScalar *xx = x->array,*y = (PetscScalar*)yin;
411:   PetscInt    i,bs = xin->bs,start,j;

413:   /*
414:        For optimization could treat bs = 2, 3, 4, 5 as special cases with loop unrolling
415:   */
417:   if (m == INSERT_VALUES) {
418:     for (i=0; i<ni; i++) {
419:       start = bs*ix[i];
420:       if (start < 0) continue;
421: #if defined(PETSC_USE_DEBUG)
422:       if (start >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",start,xin->n);
423: #endif
424:       for (j=0; j<bs; j++) {
425:         xx[start+j] = y[j];
426:       }
427:       y += bs;
428:     }
429:   } else {
430:     for (i=0; i<ni; i++) {
431:       start = bs*ix[i];
432:       if (start < 0) continue;
433: #if defined(PETSC_USE_DEBUG)
434:       if (start >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",start,xin->n);
435: #endif
436:       for (j=0; j<bs; j++) {
437:         xx[start+j] += y[j];
438:       }
439:       y += bs;
440:     }
441:   }
442:   return(0);
443: }


448: PetscErrorCode VecDestroy_Seq(Vec v)
449: {
450:   Vec_Seq        *vs = (Vec_Seq*)v->data;


455:   /* if memory was published with AMS then destroy it */
456:   PetscObjectDepublish(v);

458: #if defined(PETSC_USE_LOG)
459:   PetscLogObjectState((PetscObject)v,"Length=%D",v->n);
460: #endif
461:   if (vs->array_allocated) {PetscFree(vs->array_allocated);}
462:   PetscFree(vs);

464:   return(0);
465: }

469: static PetscErrorCode VecPublish_Seq(PetscObject obj)
470: {
472:   return(0);
473: }

475: EXTERN PetscErrorCode VecLoad_Binary(PetscViewer, VecType, Vec*);

477: static struct _VecOps DvOps = {VecDuplicate_Seq,
478:             VecDuplicateVecs_Default,
479:             VecDestroyVecs_Default,
480:             VecDot_Seq,
481:             VecMDot_Seq,
482:             VecNorm_Seq,
483:             VecTDot_Seq,
484:             VecMTDot_Seq,
485:             VecScale_Seq,
486:             VecCopy_Seq,
487:             VecSet_Seq,
488:             VecSwap_Seq,
489:             VecAXPY_Seq,
490:             VecAXPBY_Seq,
491:             VecMAXPY_Seq,
492:             VecAYPX_Seq,
493:             VecWAXPY_Seq,
494:             VecPointwiseMult_Seq,
495:             VecPointwiseDivide_Seq,
496:             VecSetValues_Seq,0,0,
497:             VecGetArray_Seq,
498:             VecGetSize_Seq,
499:             VecGetSize_Seq,
500:             VecRestoreArray_Seq,
501:             VecMax_Seq,
502:             VecMin_Seq,
503:             VecSetRandom_Seq,0,
504:             VecSetValuesBlocked_Seq,
505:             VecDestroy_Seq,
506:             VecView_Seq,
507:             VecPlaceArray_Seq,
508:             VecReplaceArray_Seq,
509:             VecDot_Seq,
510:             VecTDot_Seq,
511:             VecNorm_Seq,
512:             VecLoadIntoVector_Default,
513:             VecReciprocal_Default,
514:             0, /* VecViewNative */
515:             VecConjugate_Seq,
516:             0,
517:             0,
518:             VecResetArray_Seq,
519:             0,
520:             VecMaxPointwiseDivide_Seq,
521:             VecLoad_Binary,
522:             VecPointwiseMax_Seq,
523:             VecPointwiseMaxAbs_Seq,
524:             VecPointwiseMin_Seq,
525:             VecGetValues_Seq
526:           };


529: /*
530:       This is called by VecCreate_Seq() (i.e. VecCreateSeq()) and VecCreateSeqWithArray()
531: */
534: static PetscErrorCode VecCreate_Seq_Private(Vec v,const PetscScalar array[])
535: {
536:   Vec_Seq        *s;

540:   PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
541:   PetscNew(Vec_Seq,&s);
542:   v->data            = (void*)s;
543:   v->bops->publish   = VecPublish_Seq;
544:   v->n               = PetscMax(v->n,v->N);
545:   v->N               = PetscMax(v->n,v->N);
546:   v->petscnative     = PETSC_TRUE;
547:   s->array           = (PetscScalar *)array;
548:   s->array_allocated = 0;
549:   if (!v->map) {
550:     PetscMapCreateMPI(v->comm,v->n,v->N,&v->map);
551:   }
552:   PetscObjectChangeTypeName((PetscObject)v,VECSEQ);
553: #if defined(PETSC_HAVE_MATLAB)
554:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscMatlabEnginePut_C","VecMatlabEnginePut_Default",VecMatlabEnginePut_Default);
555:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscMatlabEngineGet_C","VecMatlabEngineGet_Default",VecMatlabEngineGet_Default);
556: #endif
557:   PetscPublishAll(v);
558:   return(0);
559: }

563: /*@C
564:    VecCreateSeqWithArray - Creates a standard,sequential array-style vector,
565:    where the user provides the array space to store the vector values.

567:    Collective on MPI_Comm

569:    Input Parameter:
570: +  comm - the communicator, should be PETSC_COMM_SELF
571: .  n - the vector length 
572: -  array - memory where the vector elements are to be stored.

574:    Output Parameter:
575: .  V - the vector

577:    Notes:
578:    Use VecDuplicate() or VecDuplicateVecs() to form additional vectors of the
579:    same type as an existing vector.

581:    If the user-provided array is PETSC_NULL, then VecPlaceArray() can be used
582:    at a later stage to SET the array for storing the vector values.

584:    PETSc does NOT free the array when the vector is destroyed via VecDestroy().
585:    The user should not free the array until the vector is destroyed.

587:    Level: intermediate

589:    Concepts: vectors^creating with array

591: .seealso: VecCreateMPIWithArray(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), 
592:           VecCreateGhost(), VecCreateSeq(), VecPlaceArray()
593: @*/
594: PetscErrorCode PETSCVEC_DLLEXPORT VecCreateSeqWithArray(MPI_Comm comm,PetscInt n,const PetscScalar array[],Vec *V)
595: {

599:   VecCreate(comm,V);
600:   VecSetSizes(*V,n,n);
601:   VecCreate_Seq_Private(*V,array);
602:   return(0);
603: }

605: /*MC
606:    VECSEQ - VECSEQ = "seq" - The basic sequential vector

608:    Options Database Keys:
609: . -vec_type seq - sets the vector type to VECSEQ during a call to VecSetFromOptions()

611:   Level: beginner

613: .seealso: VecCreate(), VecSetType(), VecSetFromOptions(), VecCreateSeqWithArray(), VECMPI, VecType, VecCreateMPI(), VecCreateSeq()
614: M*/

619: PetscErrorCode PETSCVEC_DLLEXPORT VecCreate_Seq(Vec V)
620: {
621:   Vec_Seq        *s;
622:   PetscScalar    *array;
624:   PetscInt       n = PetscMax(V->n,V->N);
625:   PetscMPIInt    size;

628:   MPI_Comm_size(V->comm,&size);
629:   if (size > 1) {
630:     SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot create VECSEQ on more than one process");
631:   }
632:   PetscMalloc( n*sizeof(PetscScalar),&array);
633:   PetscMemzero(array,n*sizeof(PetscScalar));
634:   VecCreate_Seq_Private(V,array);
635:   s    = (Vec_Seq*)V->data;
636:   s->array_allocated = array;
637:   return(0);
638: }


644: PetscErrorCode VecDuplicate_Seq(Vec win,Vec *V)
645: {

649:   VecCreateSeq(win->comm,win->n,V);
650:   if (win->mapping) {
651:     (*V)->mapping = win->mapping;
652:     PetscObjectReference((PetscObject)win->mapping);
653:   }
654:   if (win->bmapping) {
655:     (*V)->bmapping = win->bmapping;
656:     PetscObjectReference((PetscObject)win->bmapping);
657:   }
658:   (*V)->bs = win->bs;
659:   PetscOListDuplicate(win->olist,&(*V)->olist);
660:   PetscFListDuplicate(win->qlist,&(*V)->qlist);
661:   return(0);
662: }