Actual source code: bvec2.c

  1: /*
  2:    Implements the sequential vectors.
  3: */

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

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

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

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

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

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

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

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

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

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

210:   return(0);
211: }

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

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

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

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

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

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

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

298: PetscErrorCode VecView_Seq(Vec xin,PetscViewer viewer)
299: {
300:   Vec_Seq     *x = (Vec_Seq *)xin->data;
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) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
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) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
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:   } else if (issocket) {
333:     PetscViewerSocketPutScalar(viewer,xin->n,1,x->array);
334:   } else if (isbinary) {
335:     VecView_Seq_Binary(xin,viewer);
336: #if defined(PETSC_HAVE_MATHEMATICA)
337:   } else if (ismathematica) {
338:     PetscViewerMathematicaPutVector(viewer,xin);
339: #endif
340: #if defined(PETSC_HAVE_PNETCDF)
341:   } else if (isnetcdf) {
342:     VecView_Seq_Netcdf(xin,viewer);
343: #endif
344: #if defined(PETSC_HAVE_MATLAB) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
345:   } else if (ismatlab) {
346:     VecView_Seq_Matlab(xin,viewer);
347: #endif
348:   } else {
349:     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported by this vector object",((PetscObject)viewer)->type_name);
350:   }
351:   return(0);
352: }

356: PetscErrorCode VecSetValues_Seq(Vec xin,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode m)
357: {
358:   Vec_Seq     *x = (Vec_Seq *)xin->data;
359:   PetscScalar *xx = x->array;
360:   PetscInt         i;

363:   if (m == INSERT_VALUES) {
364:     for (i=0; i<ni; i++) {
365:       if (ix[i] < 0) continue;
366: #if defined(PETSC_USE_BOPT_g)
367:       if (ix[i] >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",ix[i],xin->n);
368: #endif
369:       xx[ix[i]] = y[i];
370:     }
371:   } else {
372:     for (i=0; i<ni; i++) {
373:       if (ix[i] < 0) continue;
374: #if defined(PETSC_USE_BOPT_g)
375:       if (ix[i] >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",ix[i],xin->n);
376: #endif
377:       xx[ix[i]] += y[i];
378:     }
379:   }
380:   return(0);
381: }

385: PetscErrorCode VecSetValuesBlocked_Seq(Vec xin,PetscInt ni,const PetscInt ix[],const PetscScalar yin[],InsertMode m)
386: {
387:   Vec_Seq     *x = (Vec_Seq *)xin->data;
388:   PetscScalar *xx = x->array,*y = (PetscScalar*)yin;
389:   PetscInt         i,bs = xin->bs,start,j;

391:   /*
392:        For optimization could treat bs = 2, 3, 4, 5 as special cases with loop unrolling
393:   */
395:   if (m == INSERT_VALUES) {
396:     for (i=0; i<ni; i++) {
397:       start = bs*ix[i];
398:       if (start < 0) continue;
399: #if defined(PETSC_USE_BOPT_g)
400:       if (start >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",start,xin->n);
401: #endif
402:       for (j=0; j<bs; j++) {
403:         xx[start+j] = y[j];
404:       }
405:       y += bs;
406:     }
407:   } else {
408:     for (i=0; i<ni; i++) {
409:       start = bs*ix[i];
410:       if (start < 0) continue;
411: #if defined(PETSC_USE_BOPT_g)
412:       if (start >= xin->n) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Out of range index value %D maximum %D",start,xin->n);
413: #endif
414:       for (j=0; j<bs; j++) {
415:         xx[start+j] += y[j];
416:       }
417:       y += bs;
418:     }
419:   }
420:   return(0);
421: }


426: PetscErrorCode VecDestroy_Seq(Vec v)
427: {
428:   Vec_Seq *vs = (Vec_Seq*)v->data;


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

436: #if defined(PETSC_USE_LOG)
437:   PetscLogObjectState((PetscObject)v,"Length=%D",v->n);
438: #endif
439:   if (vs->array_allocated) {PetscFree(vs->array_allocated);}
440:   PetscFree(vs);

442:   return(0);
443: }

447: static PetscErrorCode VecPublish_Seq(PetscObject obj)
448: {
450:   return(0);
451: }

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

455: static struct _VecOps DvOps = {VecDuplicate_Seq,
456:             VecDuplicateVecs_Default,
457:             VecDestroyVecs_Default,
458:             VecDot_Seq,
459:             VecMDot_Seq,
460:             VecNorm_Seq,
461:             VecTDot_Seq,
462:             VecMTDot_Seq,
463:             VecScale_Seq,
464:             VecCopy_Seq,
465:             VecSet_Seq,
466:             VecSwap_Seq,
467:             VecAXPY_Seq,
468:             VecAXPBY_Seq,
469:             VecMAXPY_Seq,
470:             VecAYPX_Seq,
471:             VecWAXPY_Seq,
472:             VecPointwiseMult_Seq,
473:             VecPointwiseDivide_Seq,
474:             VecSetValues_Seq,0,0,
475:             VecGetArray_Seq,
476:             VecGetSize_Seq,
477:             VecGetSize_Seq,
478:             VecRestoreArray_Seq,
479:             VecMax_Seq,
480:             VecMin_Seq,
481:             VecSetRandom_Seq,0,
482:             VecSetValuesBlocked_Seq,
483:             VecDestroy_Seq,
484:             VecView_Seq,
485:             VecPlaceArray_Seq,
486:             VecReplaceArray_Seq,
487:             VecDot_Seq,
488:             VecTDot_Seq,
489:             VecNorm_Seq,
490:             VecLoadIntoVector_Default,
491:             VecReciprocal_Default,
492:             0, /* VecViewNative */
493:             VecConjugate_Seq,
494:             0,
495:             0,
496:             VecResetArray_Seq,
497:             0,
498:             VecMaxPointwiseDivide_Seq,
499:             VecLoad_Binary};


502: /*
503:       This is called by VecCreate_Seq() (i.e. VecCreateSeq()) and VecCreateSeqWithArray()
504: */
507: static PetscErrorCode VecCreate_Seq_Private(Vec v,const PetscScalar array[])
508: {
509:   Vec_Seq *s;

513:   PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
514:   PetscNew(Vec_Seq,&s);
515:   PetscMemzero(s,sizeof(Vec_Seq));
516:   v->data            = (void*)s;
517:   v->bops->publish   = VecPublish_Seq;
518:   v->n               = PetscMax(v->n,v->N);
519:   v->N               = PetscMax(v->n,v->N);
520:   v->petscnative     = PETSC_TRUE;
521:   s->array           = (PetscScalar *)array;
522:   s->array_allocated = 0;
523:   if (!v->map) {
524:     PetscMapCreateMPI(v->comm,v->n,v->N,&v->map);
525:   }
526:   PetscObjectChangeTypeName((PetscObject)v,VECSEQ);
527: #if defined(PETSC_HAVE_MATLAB) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
528:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscMatlabEnginePut_C","VecMatlabEnginePut_Default",VecMatlabEnginePut_Default);
529:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscMatlabEngineGet_C","VecMatlabEngineGet_Default",VecMatlabEngineGet_Default);
530: #endif
531:   PetscPublishAll(v);
532:   return(0);
533: }

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

541:    Collective on MPI_Comm

543:    Input Parameter:
544: +  comm - the communicator, should be PETSC_COMM_SELF
545: .  n - the vector length 
546: -  array - memory where the vector elements are to be stored.

548:    Output Parameter:
549: .  V - the vector

551:    Notes:
552:    Use VecDuplicate() or VecDuplicateVecs() to form additional vectors of the
553:    same type as an existing vector.

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

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

561:    Level: intermediate

563:    Concepts: vectors^creating with array

565: .seealso: VecCreateMPIWithArray(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), 
566:           VecCreateGhost(), VecCreateSeq(), VecPlaceArray()
567: @*/
568: PetscErrorCode VecCreateSeqWithArray(MPI_Comm comm,PetscInt n,const PetscScalar array[],Vec *V)
569: {

573:   VecCreate(comm,V);
574:   VecSetSizes(*V,n,n);
575:   VecCreate_Seq_Private(*V,array);
576:   return(0);
577: }

579: /*MC
580:    VECSEQ - VECSEQ = "seq" - The basic sequential vector

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

585:   Level: beginner

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

593: PetscErrorCode VecCreate_Seq(Vec V)
594: {
595:   Vec_Seq        *s;
596:   PetscScalar    *array;
598:   PetscInt            n = PetscMax(V->n,V->N);
599:   PetscMPIInt    size;

602:   MPI_Comm_size(V->comm,&size);
603:   if (size > 1) {
604:     SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot create VECSEQ on more than one process");
605:   }
606:   PetscMalloc( n*sizeof(PetscScalar),&array);
607:   PetscMemzero(array,n*sizeof(PetscScalar));
608:   VecCreate_Seq_Private(V,array);
609:   s    = (Vec_Seq*)V->data;
610:   s->array_allocated = array;
611:   return(0);
612: }


618: PetscErrorCode VecDuplicate_Seq(Vec win,Vec *V)
619: {

623:   VecCreateSeq(win->comm,win->n,V);
624:   if (win->mapping) {
625:     (*V)->mapping = win->mapping;
626:     PetscObjectReference((PetscObject)win->mapping);
627:   }
628:   if (win->bmapping) {
629:     (*V)->bmapping = win->bmapping;
630:     PetscObjectReference((PetscObject)win->bmapping);
631:   }
632:   (*V)->bs = win->bs;
633:   PetscOListDuplicate(win->olist,&(*V)->olist);
634:   PetscFListDuplicate(win->qlist,&(*V)->qlist);
635:   return(0);
636: }