Actual source code: binv.c

  1: #define PETSC_DLL
 2:  #include src/sys/viewer/viewerimpl.h
 3:  #include petscsys.h
  4: #include <fcntl.h>
  5: #if defined(PETSC_HAVE_UNISTD_H)
  6: #include <unistd.h>
  7: #endif
  8: #if defined (PETSC_HAVE_IO_H)
  9: #include <io.h>
 10: #endif

 12: typedef struct  {
 13:   int           fdes;            /* file descriptor */
 14:   PetscFileMode btype;           /* read or write? */
 15:   FILE          *fdes_info;      /* optional file containing info on binary file*/
 16:   PetscTruth    storecompressed; /* gzip the write binary file when closing it*/
 17:   char          *filename;
 18:   PetscTruth    skipinfo;        /* Don't create info file for writing; don't use for reading */
 19:   PetscTruth    skipoptions;     /* don't use PETSc options database when loading */
 20: } PetscViewer_Binary;

 24: PetscErrorCode PetscViewerGetSingleton_Binary(PetscViewer viewer,PetscViewer *outviewer)
 25: {
 26:   int                rank;
 27:   PetscErrorCode     ierr;
 28:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data,*obinary;

 31:   MPI_Comm_rank(viewer->comm,&rank);
 32:   if (!rank) {
 33:     PetscViewerCreate(PETSC_COMM_SELF,outviewer);
 34:     PetscViewerSetType(*outviewer,PETSC_VIEWER_BINARY);
 35:     obinary = (PetscViewer_Binary*)(*outviewer)->data;
 36:     PetscMemcpy(obinary,vbinary,sizeof(PetscViewer_Binary));
 37:   } else {
 38:     *outviewer = 0;
 39:   }
 40:   return(0);
 41: }

 45: PetscErrorCode PetscViewerRestoreSingleton_Binary(PetscViewer viewer,PetscViewer *outviewer)
 46: {
 48:   PetscErrorCode rank;

 51:   MPI_Comm_rank(viewer->comm,&rank);
 52:   if (!rank) {
 53:     PetscFree((*outviewer)->data);
 54:     PetscHeaderDestroy(*outviewer);
 55:   }
 56:   return(0);
 57: }

 61: /*@C
 62:     PetscViewerBinaryGetDescriptor - Extracts the file descriptor from a PetscViewer.

 64:     Not Collective

 66: +   viewer - PetscViewer context, obtained from PetscViewerBinaryOpen()
 67: -   fdes - file descriptor

 69:     Level: advanced

 71:     Notes:
 72:       For writable binary PetscViewers, the descriptor will only be valid for the 
 73:     first processor in the communicator that shares the PetscViewer. For readable 
 74:     files it will only be valid on nodes that have the file. If node 0 does not
 75:     have the file it generates an error even if another node does have the file.
 76:  
 77:     Fortran Note:
 78:     This routine is not supported in Fortran.

 80:   Concepts: file descriptor^getting
 81:   Concepts: PetscViewerBinary^accessing file descriptor

 83: .seealso: PetscViewerBinaryOpen(),PetscViewerBinaryGetInfoPointer()
 84: @*/
 85: PetscErrorCode  PetscViewerBinaryGetDescriptor(PetscViewer viewer,int *fdes)
 86: {
 87:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

 90:   *fdes = vbinary->fdes;
 91:   return(0);
 92: }

 96: /*@
 97:     PetscViewerBinarySkipInfo - Binary file will not have .info file created with it

 99:     Not Collective

101:     Input Paramter:
102: .   viewer - PetscViewer context, obtained from PetscViewerBinaryOpen()

104:     Options Database Key:
105: .   -viewer_binary_skip_info

107:     Level: advanced

109:     Notes: This must be called after PetscViewerSetType() but before PetscViewerBinarySetFilename()

111:    Concepts: PetscViewerBinary^accessing info file

113: .seealso: PetscViewerBinaryOpen(), PetscViewerBinaryGetDescriptor(), PetscViewerBinarySetSkipOptions(),
114:           PetscViewerBinaryGetSkipOptions()
115: @*/
116: PetscErrorCode  PetscViewerBinarySkipInfo(PetscViewer viewer)
117: {
118:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

121:   vbinary->skipinfo = PETSC_TRUE;
122:   return(0);
123: }

127: /*@
128:     PetscViewerBinarySetSkipOptions - do not use the PETSc options database when loading objects

130:     Not Collective

132:     Input Paramter:
133: +   viewer - PetscViewer context, obtained from PetscViewerBinaryOpen()
134: -   skip - PETSC_TRUE means do not use

136:     Options Database Key:
137: .   -viewer_binary_skip_options

139:     Level: advanced

141:     Notes: This must be called after PetscViewerSetType()

143:    Concepts: PetscViewerBinary^accessing info file

145: .seealso: PetscViewerBinaryOpen(), PetscViewerBinaryGetDescriptor(), PetscViewerBinarySkipInfo(),
146:           PetscViewerBinaryGetSkipOptions()
147: @*/
148: PetscErrorCode  PetscViewerBinarySetSkipOptions(PetscViewer viewer,PetscTruth skip)
149: {
150:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

153:   vbinary->skipoptions = skip;
154:   return(0);
155: }

159: /*@
160:     PetscViewerBinaryGetSkipOptions - checks if viewer uses the PETSc options database when loading objects

162:     Not Collective

164:     Input Paramter:
165: .   viewer - PetscViewer context, obtained from PetscViewerBinaryOpen()

167:     Output Parameter:
168: .   skip - PETSC_TRUE means do not use

170:     Level: advanced

172:     Notes: This must be called after PetscViewerSetType()

174:    Concepts: PetscViewerBinary^accessing info file

176: .seealso: PetscViewerBinaryOpen(), PetscViewerBinaryGetDescriptor(), PetscViewerBinarySkipInfo(),
177:           PetscViewerBinarySetSkipOptions()
178: @*/
179: PetscErrorCode  PetscViewerBinaryGetSkipOptions(PetscViewer viewer,PetscTruth *skip)
180: {
181:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

184:   *skip = vbinary->skipoptions;
185:   return(0);
186: }

190: /*@C
191:     PetscViewerBinaryGetInfoPointer - Extracts the file pointer for the ASCII
192:           info file associated with a binary file.

194:     Not Collective

196: +   viewer - PetscViewer context, obtained from PetscViewerBinaryOpen()
197: -   file - file pointer

199:     Level: advanced

201:     Notes:
202:       For writable binary PetscViewers, the descriptor will only be valid for the 
203:     first processor in the communicator that shares the PetscViewer.
204:  
205:     Fortran Note:
206:     This routine is not supported in Fortran.

208:   Concepts: PetscViewerBinary^accessing info file

210: .seealso: PetscViewerBinaryOpen(),PetscViewerBinaryGetDescriptor()
211: @*/
212: PetscErrorCode  PetscViewerBinaryGetInfoPointer(PetscViewer viewer,FILE **file)
213: {
214:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

217:   *file = vbinary->fdes_info;
218:   return(0);
219: }

223: PetscErrorCode PetscViewerDestroy_Binary(PetscViewer v)
224: {
225:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)v->data;
226:   PetscErrorCode     ierr;
227:   int                rank;

230:   MPI_Comm_rank(v->comm,&rank);
231:   if ((!rank || vbinary->btype == FILE_MODE_READ) && vbinary->fdes) {
232:     close(vbinary->fdes);
233:     if (!rank && vbinary->storecompressed) {
234:       char par[PETSC_MAX_PATH_LEN],buf[PETSC_MAX_PATH_LEN];
235:       FILE *fp;
236:       /* compress the file */
237:       PetscStrcpy(par,"gzip ");
238:       PetscStrcat(par,vbinary->filename);
239: #if defined(PETSC_HAVE_POPEN)
240:       PetscPOpen(PETSC_COMM_SELF,PETSC_NULL,par,"r",&fp);
241:       if (fgets(buf,1024,fp)) {
242:         SETERRQ2(PETSC_ERR_LIB,"Error from command %s\n%s",par,buf);
243:       }
244:       PetscPClose(PETSC_COMM_SELF,fp);
245: #else
246:       SETERRQ(PETSC_ERR_SUP_SYS,"Cannot run external programs on this machine");
247: #endif
248:     }
249:   }
250:   if (vbinary->fdes_info) fclose(vbinary->fdes_info);
251:   PetscStrfree(vbinary->filename);
252:   PetscFree(vbinary);
253:   return(0);
254: }

258: /*@
259:    PetscViewerBinaryCreate - Create a binary viewer.

261:    Collective on MPI_Comm

263:    Input Parameters:
264: .  comm - MPI communicator

266:    Output Parameter:
267: .  binv - PetscViewer for binary input/output

269:    Level: beginner
270: @*/
271: PetscErrorCode  PetscViewerBinaryCreate(MPI_Comm comm,PetscViewer *binv)
272: {
274: 
276:   PetscViewerCreate(comm,binv);
277:   PetscViewerSetType(*binv,PETSC_VIEWER_BINARY);
278:   return(0);
279: }

283: /*@C
284:    PetscViewerBinaryOpen - Opens a file for binary input/output.

286:    Collective on MPI_Comm

288:    Input Parameters:
289: +  comm - MPI communicator
290: .  name - name of file 
291: -  type - type of file
292: $    FILE_MODE_WRITE - create new file for binary output
293: $    FILE_MODE_READ - open existing file for binary input
294: $    FILE_MODE_APPEND - open existing file for binary output

296:    Output Parameter:
297: .  binv - PetscViewer for binary input/output to use with the specified file

299:    Level: beginner

301:    Note:
302:    This PetscViewer should be destroyed with PetscViewerDestroy().

304:     For reading files, the filename may begin with ftp:// or http:// and/or
305:     end with .gz; in this case file is brought over and uncompressed.

307:     For creating files, if the file name ends with .gz it is automatically 
308:     compressed when closed.

310:     For writing files it only opens the file on processor 0 in the communicator.
311:     For readable files it opens the file on all nodes that have the file. If 
312:     node 0 does not have the file it generates an error even if other nodes
313:     do have the file.

315:    Concepts: binary files
316:    Concepts: PetscViewerBinary^creating
317:    Concepts: gzip
318:    Concepts: accessing remote file
319:    Concepts: remote file

321: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
322:           VecView(), MatView(), VecLoad(), MatLoad(), PetscViewerBinaryGetDescriptor(),
323:           PetscViewerBinaryGetInfoPointer(), PetscFileMode, PetscViewer, PetscViewerBinaryRead()
324: @*/
325: PetscErrorCode  PetscViewerBinaryOpen(MPI_Comm comm,const char name[],PetscFileMode type,PetscViewer *binv)
326: {
328: 
330:   PetscViewerCreate(comm,binv);
331:   PetscViewerSetType(*binv,PETSC_VIEWER_BINARY);
332:   PetscViewerFileSetMode(*binv,type);
333:   PetscViewerFileSetName(*binv,name);
334:   return(0);
335: }

339: /*@C
340:    PetscViewerBinaryRead - Reads from a binary file, all processors get the same result

342:    Collective on MPI_Comm

344:    Input Parameters:
345: +  viewer - the binary viewer
346: .  data - location to write the data
347: .  count - number of items of data to read
348: -  datatype - type of data to read

350:    Level: beginner

352:    Concepts: binary files

354: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
355:           VecView(), MatView(), VecLoad(), MatLoad(), PetscViewerBinaryGetDescriptor(),
356:           PetscViewerBinaryGetInfoPointer(), PetscFileMode, PetscViewer, PetscBinaryViewerRead()
357: @*/
358: PetscErrorCode  PetscViewerBinaryRead(PetscViewer viewer,void *data,PetscInt count,PetscDataType dtype)
359: {
360:   PetscErrorCode     ierr;
361:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

363:   PetscSynchronizedBinaryRead(viewer->comm,vbinary->fdes,data,count,dtype);
364:   return(0);
365: }

369: /*@C
370:    PetscViewerBinaryWrite - writes to a binary file, only from the first process

372:    Collective on MPI_Comm

374:    Input Parameters:
375: +  viewer - the binary viewer
376: .  data - location of data
377: .  count - number of items of data to read
378: .  istemp - data may be overwritten
379: -  datatype - type of data to read

381:    Level: beginner

383:    Concepts: binary files

385: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
386:           VecView(), MatView(), VecLoad(), MatLoad(), PetscViewerBinaryGetDescriptor(),
387:           PetscViewerBinaryGetInfoPointer(), PetscFileMode, PetscViewer, PetscBinaryViewerRead()
388: @*/
389: PetscErrorCode  PetscViewerBinaryWrite(PetscViewer viewer,void *data,PetscInt count,PetscDataType dtype,PetscTruth istemp)
390: {
391:   PetscErrorCode     ierr;
392:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

394:   PetscSynchronizedBinaryWrite(viewer->comm,vbinary->fdes,data,count,dtype,istemp);
395:   return(0);
396: }

400: /*@C
401:    PetscViewerBinaryWriteStringArray - writes to a binary file, only from the first process an array of strings

403:    Collective on MPI_Comm

405:    Input Parameters:
406: +  viewer - the binary viewer
407: -  data - location of the array of strings


410:    Level: intermediate

412:    Concepts: binary files

414:     Notes: array of strings is null terminated

416: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
417:           VecView(), MatView(), VecLoad(), MatLoad(), PetscViewerBinaryGetDescriptor(),
418:           PetscViewerBinaryGetInfoPointer(), PetscFileMode, PetscViewer, PetscBinaryViewerRead()
419: @*/
420: PetscErrorCode  PetscViewerBinaryWriteStringArray(PetscViewer viewer,char **data)
421: {
422:   PetscErrorCode     ierr;
423:   PetscInt           i,n = 0,*sizes;

425:   /* count number of strings */
426:   while (data[n++]);
427:   n--;
428:   PetscMalloc((n+1)*sizeof(PetscInt),&sizes);
429:   sizes[0] = n;
430:   for (i=0; i<n; i++) {
431:     size_t tmp;
432:     PetscStrlen(data[i],&tmp);
433:     sizes[i+1] = tmp + 1;   /* size includes space for the null terminator */
434:   }
435:   PetscViewerBinaryWrite(viewer,sizes,n+1,PETSC_INT,PETSC_FALSE);
436:   for (i=0; i<n; i++) {
437:     PetscViewerBinaryWrite(viewer,data[i],sizes[i+1],PETSC_CHAR,PETSC_FALSE);
438:   }
439:   PetscFree(sizes);
440:   return(0);
441: }

443: /*@C
444:    PetscViewerBinaryReadStringArray - reads a binary file an array of strings

446:    Collective on MPI_Comm

448:    Input Parameter:
449: .  viewer - the binary viewer

451:    Output Parameter:
452: .  data - location of the array of strings

454:    Level: intermediate

456:    Concepts: binary files

458:     Notes: array of strings is null terminated

460: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
461:           VecView(), MatView(), VecLoad(), MatLoad(), PetscViewerBinaryGetDescriptor(),
462:           PetscViewerBinaryGetInfoPointer(), PetscFileMode, PetscViewer, PetscBinaryViewerRead()
463: @*/
464: PetscErrorCode  PetscViewerBinaryReadStringArray(PetscViewer viewer,char ***data)
465: {
466:   PetscErrorCode     ierr;
467:   PetscInt           i,n,*sizes,N = 0;

469:   /* count number of strings */
470:   PetscViewerBinaryRead(viewer,&n,1,PETSC_INT);
471:   PetscMalloc(n*sizeof(PetscInt),&sizes);
472:   PetscViewerBinaryRead(viewer,sizes,n,PETSC_INT);
473:   for (i=0; i<n; i++) {
474:     N += sizes[i];
475:   }
476:   PetscMalloc((n+1)*sizeof(char*) + N*sizeof(char),data);
477:   (*data)[0] = (char*)((*data) + n + 1);
478:   for (i=1; i<n; i++) {
479:     (*data)[i] = (*data)[i-1] + sizes[i-1];
480:   }
481:   PetscViewerBinaryRead(viewer,(*data)[0],N,PETSC_CHAR);
482:   (*data)[n] = 0;
483:   PetscFree(sizes);
484:   return(0);
485: }

489: /*@C
490:      PetscViewerFileGetMode - Gets the type of file to be open

492:     Collective on PetscViewer

494:   Input Parameter:
495: .  viewer - the PetscViewer; must be a binary, Matlab, hdf, or netcdf PetscViewer

497:   Output Parameter:
498: .  type - type of file
499: $    FILE_MODE_WRITE - create new file for binary output
500: $    FILE_MODE_READ - open existing file for binary input
501: $    FILE_MODE_APPEND - open existing file for binary output

503:   Level: advanced

505: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen()

507: @*/
508: PetscErrorCode  PetscViewerFileGetMode(PetscViewer viewer,PetscFileMode *type)
509: {
510:   PetscErrorCode ierr,(*f)(PetscViewer,PetscFileMode*);

515:   PetscObjectQueryFunction((PetscObject)viewer,"PetscViewerFileGetMode_C",(void (**)(void))&f);
516:   if (f) {
517:     (*f)(viewer,type);
518:   }
519:   return(0);
520: }

524: /*@C
525:      PetscViewerFileSetMode - Sets the type of file to be open

527:     Collective on PetscViewer

529:   Input Parameters:
530: +  viewer - the PetscViewer; must be a binary, Matlab, hdf, or netcdf PetscViewer
531: -  type - type of file
532: $    FILE_MODE_WRITE - create new file for binary output
533: $    FILE_MODE_READ - open existing file for binary input
534: $    FILE_MODE_APPEND - open existing file for binary output

536:   Level: advanced

538: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen()

540: @*/
541: PetscErrorCode  PetscViewerFileSetMode(PetscViewer viewer,PetscFileMode type)
542: {
543:   PetscErrorCode ierr,(*f)(PetscViewer,PetscFileMode);

547:   PetscObjectQueryFunction((PetscObject)viewer,"PetscViewerFileSetMode_C",(void (**)(void))&f);
548:   if (f) {
549:     (*f)(viewer,type);
550:   }
551:   return(0);
552: }

557: PetscErrorCode  PetscViewerFileGetMode_Binary(PetscViewer viewer,PetscFileMode *type)
558: {
559:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

562:   *type = vbinary->btype;
563:   return(0);
564: }

570: PetscErrorCode  PetscViewerFileSetMode_Binary(PetscViewer viewer,PetscFileMode type)
571: {
572:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

575:   vbinary->btype = type;
576:   return(0);
577: }

582: /*@
583:     PetscViewerBinaryLoadInfo - Loads options from the name.info file
584:        if it exists.

586:    Collective on PetscViewer

588:   Input Parameter:
589: .    viewer - the binary viewer whose options you wish to load

591:    Level: developer

593: @*/
594: PetscErrorCode  PetscViewerBinaryLoadInfo(PetscViewer viewer)
595: {
596:   FILE               *file;
597:   char               string[256],*first,*second,*final;
598:   size_t             len;
599:   PetscErrorCode     ierr;
600:   PetscToken         *token;
601:   PetscViewer_Binary *vbinary = (PetscViewer_Binary*)viewer->data;

604:   if (vbinary->skipinfo) return(0);

606:   PetscViewerBinaryGetInfoPointer(viewer,&file);
607:   if (!file) return(0);

609:   /* read rows of the file adding them to options database */
610:   while (fgets(string,256,file)) {
611:     /* Comments are indicated by #, ! or % in the first column */
612:     if (string[0] == '#') continue;
613:     if (string[0] == '!') continue;
614:     if (string[0] == '%') continue;
615:     PetscTokenCreate(string,' ',&token);
616:     PetscTokenFind(token,&first);
617:     PetscTokenFind(token,&second);
618:     if (first && first[0] == '-') {
619:       PetscTruth wrongtype;
620:       /*
621:          Check for -mat_complex or -mat_double
622:       */
623: #if defined(PETSC_USE_COMPLEX)
624:       PetscStrncmp(first,"-mat_double",11,&wrongtype);
625:       if (wrongtype) {
626:         SETERRQ(PETSC_ERR_FILE_UNEXPECTED,"Loading double number matrix with complex number code");
627:       }
628: #else
629:       PetscStrncmp(first,"-mat_complex",12,&wrongtype);
630:       if (wrongtype) {
631:         SETERRQ(PETSC_ERR_FILE_UNEXPECTED,"Loading complex number matrix with double number code");
632:       }
633: #endif

635:       if (second) {final = second;} else {final = first;}
636:       PetscStrlen(final,&len);
637:       while (len > 0 && (final[len-1] == ' ' || final[len-1] == '\n')) {
638:         len--; final[len] = 0;
639:       }
640:       PetscOptionsSetValue(first,second);
641:     }
642:     PetscTokenDestroy(token);
643:   }
644:   return(0);
645: }

647: /*
648:         Actually opens the file 
649: */
653: PetscErrorCode  PetscViewerFileSetName_Binary(PetscViewer viewer,const char name[])
654: {
655:   int                 rank;
656:   PetscErrorCode      ierr;
657:   size_t              len;
658:   PetscViewer_Binary  *vbinary = (PetscViewer_Binary*)viewer->data;
659:   const char          *fname;
660:   char                bname[PETSC_MAX_PATH_LEN],*gz;
661:   PetscTruth          found;
662:   PetscFileMode type = vbinary->btype;

665:   if (type == (PetscFileMode) -1) {
666:     SETERRQ(PETSC_ERR_ORDER,"Must call PetscViewerBinarySetFileType() before PetscViewerFileSetName()");
667:   }
668:   PetscOptionsGetTruth(viewer->prefix,"-viewer_binary_skip_info",&vbinary->skipinfo,PETSC_NULL);
669:   PetscOptionsGetTruth(viewer->prefix,"-viewer_binary_skip_options",&vbinary->skipoptions,PETSC_NULL);

671:   MPI_Comm_rank(viewer->comm,&rank);

673:   /* copy name so we can edit it */
674:   PetscStrallocpy(name,&vbinary->filename);

676:   /* if ends in .gz strip that off and note user wants file compressed */
677:   vbinary->storecompressed = PETSC_FALSE;
678:   if (!rank && type == FILE_MODE_WRITE) {
679:     /* remove .gz if it ends library name */
680:     PetscStrstr(vbinary->filename,".gz",&gz);
681:     if (gz) {
682:       PetscStrlen(gz,&len);
683:       if (len == 3) {
684:         *gz = 0;
685:         vbinary->storecompressed = PETSC_TRUE;
686:       }
687:     }
688:   }

690:   /* only first processor opens file if writeable */
691:   if (!rank || type == FILE_MODE_READ) {

693:     if (type == FILE_MODE_READ){
694:       /* possibly get the file from remote site or compressed file */
695:       PetscFileRetrieve(viewer->comm,vbinary->filename,bname,PETSC_MAX_PATH_LEN,&found);
696:       fname = bname;
697:       if (!rank && !found) {
698:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot locate file: %s on node zero",vbinary->filename);
699:       } else if (!found) {
700:         PetscInfo(viewer,"Nonzero processor did not locate readonly file\n");
701:         fname = 0;
702:       }
703:     } else {
704:       fname = vbinary->filename;
705:     }

707: #if defined(PETSC_HAVE_O_BINARY)
708:     if (type == FILE_MODE_WRITE) {
709:       if ((vbinary->fdes = open(fname,O_WRONLY|O_CREAT|O_TRUNC|O_BINARY,0666)) == -1) {
710:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot create file %s for writing",fname);
711:       }
712:     } else if (type == FILE_MODE_READ && fname) {
713:       if ((vbinary->fdes = open(fname,O_RDONLY|O_BINARY,0)) == -1) {
714:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open file %s for reading",fname);
715:       }
716:     } else if (type == FILE_MODE_APPEND) {
717:       if ((vbinary->fdes = open(fname,O_WRONLY|O_BINARY,0)) == -1) {
718:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open file %s for writing",fname);
719:       }
720:     } else if (fname) {
721:       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown file type");
722:     }
723: #else
724:     if (type == FILE_MODE_WRITE) {
725:       if ((vbinary->fdes = creat(fname,0666)) == -1) {
726:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot create file %s for writing",fname);
727:       }
728:     } else if (type == FILE_MODE_READ && fname) {
729:       if ((vbinary->fdes = open(fname,O_RDONLY,0)) == -1) {
730:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open file %s for reading",fname);
731:       }
732:     } else if (type == FILE_MODE_APPEND) {
733:       if ((vbinary->fdes = open(fname,O_WRONLY|O_APPEND,0)) == -1) {
734:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open file %s for writing",fname);
735:       }
736:     } else if (fname) {
737:       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown file type");
738:     }
739: #endif
740:   } else vbinary->fdes = -1;
741:   viewer->format = PETSC_VIEWER_NOFORMAT;

743:   /* 
744:       try to open info file: all processors open this file if read only
745:   */
746:   if (!rank || type == FILE_MODE_READ) {
747:     char infoname[PETSC_MAX_PATH_LEN],iname[PETSC_MAX_PATH_LEN];
748: 
749:     PetscStrcpy(infoname,name);
750:     /* remove .gz if it ends library name */
751:     PetscStrstr(infoname,".gz",&gz);
752:     if (gz) {
753:       PetscStrlen(gz,&len);
754:       if (len == 3) {
755:         *gz = 0;
756:       }
757:     }
758: 
759:     PetscStrcat(infoname,".info");
760:     PetscFixFilename(infoname,iname);
761:     if (type == FILE_MODE_READ) {
762:       PetscFileRetrieve(viewer->comm,iname,infoname,PETSC_MAX_PATH_LEN,&found);
763:       if (found) {
764:         vbinary->fdes_info = fopen(infoname,"r");
765:         if (vbinary->fdes_info) {
766:           PetscViewerBinaryLoadInfo(viewer);
767:           fclose(vbinary->fdes_info);
768:         }
769:         vbinary->fdes_info = fopen(infoname,"r");
770:       }
771:     } else if (!vbinary->skipinfo) {
772:       vbinary->fdes_info = fopen(infoname,"w");
773:       if (!vbinary->fdes_info) {
774:         SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open .info file %s for writing",infoname);
775:       }
776:     }
777:   }

779: #if defined(PETSC_USE_LOG)
780:   PetscLogObjectState((PetscObject)viewer,"File: %s",name);
781: #endif
782:   return(0);
783: }

789: PetscErrorCode  PetscViewerCreate_Binary(PetscViewer v)
790: {
791:   PetscErrorCode     ierr;
792:   PetscViewer_Binary *vbinary;

795:   PetscNew(PetscViewer_Binary,&vbinary);
796:   v->data            = (void*)vbinary;
797:   v->ops->destroy    = PetscViewerDestroy_Binary;
798:   v->ops->flush      = 0;
799:   v->iformat         = 0;
800:   vbinary->fdes_info = 0;
801:   vbinary->fdes      = 0;
802:   vbinary->skipinfo        = PETSC_FALSE;
803:   vbinary->skipoptions     = PETSC_TRUE;
804:   v->ops->getsingleton     = PetscViewerGetSingleton_Binary;
805:   v->ops->restoresingleton = PetscViewerRestoreSingleton_Binary;
806:   vbinary->btype           = (PetscFileMode) -1;
807:   vbinary->storecompressed = PETSC_FALSE;
808:   vbinary->filename        = 0;

810:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetName_C",
811:                                     "PetscViewerFileSetName_Binary",
812:                                      PetscViewerFileSetName_Binary);
813:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetMode_C",
814:                                     "PetscViewerFileSetMode_Binary",
815:                                      PetscViewerFileSetMode_Binary);
816:   PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileGetMode_C",
817:                                     "PetscViewerFileGetMode_Binary",
818:                                      PetscViewerFileGetMode_Binary);
819:   return(0);
820: }


824: /* ---------------------------------------------------------------------*/
825: /*
826:     The variable Petsc_Viewer_Binary_keyval is used to indicate an MPI attribute that
827:   is attached to a communicator, in this case the attribute is a PetscViewer.
828: */
829: static int Petsc_Viewer_Binary_keyval = MPI_KEYVAL_INVALID;

833: /*@C
834:      PETSC_VIEWER_BINARY_ - Creates a binary PetscViewer shared by all processors 
835:                      in a communicator.

837:      Collective on MPI_Comm

839:      Input Parameter:
840: .    comm - the MPI communicator to share the binary PetscViewer
841:     
842:      Level: intermediate

844:    Options Database Keys:
845: $    -viewer_binary_filename <name>

847:    Environmental variables:
848: -   PETSC_VIEWER_BINARY_FILENAME

850:      Notes:
851:      Unlike almost all other PETSc routines, PETSC_VIEWER_BINARY_ does not return 
852:      an error code.  The binary PetscViewer is usually used in the form
853: $       XXXView(XXX object,PETSC_VIEWER_BINARY_(comm));

855: .seealso: PETSC_VIEWER_BINARY_WORLD, PETSC_VIEWER_BINARY_SELF, PetscViewerBinaryOpen(), PetscViewerCreate(),
856:           PetscViewerDestroy()
857: @*/
858: PetscViewer  PETSC_VIEWER_BINARY_(MPI_Comm comm)
859: {
861:   PetscTruth     flg;
862:   PetscViewer    viewer;
863:   char           fname[PETSC_MAX_PATH_LEN];

866:   if (Petsc_Viewer_Binary_keyval == MPI_KEYVAL_INVALID) {
867:     MPI_Keyval_create(MPI_NULL_COPY_FN,MPI_NULL_DELETE_FN,&Petsc_Viewer_Binary_keyval,0);
868:     if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
869:   }
870:   MPI_Attr_get(comm,Petsc_Viewer_Binary_keyval,(void **)&viewer,(int*)&flg);
871:   if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
872:   if (!flg) { /* PetscViewer not yet created */
873:     PetscOptionsGetenv(comm,"PETSC_VIEWER_BINARY_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
874:     if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
875:     if (!flg) {
876:       PetscStrcpy(fname,"binaryoutput");
877:       if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
878:     }
879:     PetscViewerBinaryOpen(comm,fname,FILE_MODE_WRITE,&viewer);
880:     if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
881:     PetscObjectRegisterDestroy((PetscObject)viewer);
882:     if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
883:     MPI_Attr_put(comm,Petsc_Viewer_Binary_keyval,(void*)viewer);
884:     if (ierr) {PetscError(__LINE__,"PETSC_VIEWER_BINARY_",__FILE__,__SDIR__,1,1," ");return(0);}
885:   }
886:   PetscFunctionReturn(viewer);
887: }