1235f7792SMatthew G. Knepley #include <petscis.h> /*I "petscis.h" I*/ 2af0996ceSBarry Smith #include <petsc/private/isimpl.h> 390e28553SVaclav Hapla #include <petsc/private/viewerimpl.h> 420e823e8SBarry Smith #include <petsclayouthdf5.h> 5235f7792SMatthew G. Knepley 6f92d9ac3SLisandro Dalcin PetscErrorCode ISView_Binary(IS is,PetscViewer viewer) 7f92d9ac3SLisandro Dalcin { 8f92d9ac3SLisandro Dalcin PetscBool skipHeader; 9f92d9ac3SLisandro Dalcin PetscLayout map; 10f92d9ac3SLisandro Dalcin PetscInt tr[2],n,s,N; 11f92d9ac3SLisandro Dalcin const PetscInt *iarray; 12f92d9ac3SLisandro Dalcin 13f92d9ac3SLisandro Dalcin PetscFunctionBegin; 14*9566063dSJacob Faibussowitsch PetscCall(PetscViewerSetUp(viewer)); 15*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryGetSkipHeader(viewer,&skipHeader)); 16f92d9ac3SLisandro Dalcin 17*9566063dSJacob Faibussowitsch PetscCall(ISGetLayout(is,&map)); 18*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetLocalSize(map,&n)); 19*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetRange(map,&s,NULL)); 20*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetSize(map,&N)); 21f92d9ac3SLisandro Dalcin 22f92d9ac3SLisandro Dalcin /* write IS header */ 23f92d9ac3SLisandro Dalcin tr[0] = IS_FILE_CLASSID; tr[1] = N; 24*9566063dSJacob Faibussowitsch if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer,tr,2,PETSC_INT)); 25f92d9ac3SLisandro Dalcin 26f92d9ac3SLisandro Dalcin /* write IS indices */ 27*9566063dSJacob Faibussowitsch PetscCall(ISGetIndices(is,&iarray)); 28*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryWriteAll(viewer,iarray,n,s,N,PETSC_INT)); 29*9566063dSJacob Faibussowitsch PetscCall(ISRestoreIndices(is,&iarray)); 30f92d9ac3SLisandro Dalcin PetscFunctionReturn(0); 31f92d9ac3SLisandro Dalcin } 32f92d9ac3SLisandro Dalcin 33235f7792SMatthew G. Knepley #if defined(PETSC_HAVE_HDF5) 34235f7792SMatthew G. Knepley /* 35235f7792SMatthew G. Knepley This should handle properly the cases where PetscInt is 32 or 64 and hsize_t is 32 or 64. These means properly casting with 36235f7792SMatthew G. Knepley checks back and forth between the two types of variables. 37235f7792SMatthew G. Knepley */ 38235f7792SMatthew G. Knepley PetscErrorCode ISLoad_HDF5(IS is, PetscViewer viewer) 39235f7792SMatthew G. Knepley { 40235f7792SMatthew G. Knepley hid_t inttype; /* int type (H5T_NATIVE_INT or H5T_NATIVE_LLONG) */ 4190e28553SVaclav Hapla PetscInt *ind; 42235f7792SMatthew G. Knepley const char *isname; 43235f7792SMatthew G. Knepley 44235f7792SMatthew G. Knepley PetscFunctionBegin; 45fa7c0e95SVaclav Hapla PetscCheck(((PetscObject)is)->name, PetscObjectComm((PetscObject)is), PETSC_ERR_SUP, "IS name must be given using PetscObjectSetName() before ISLoad() since HDF5 can store multiple objects in a single file"); 46235f7792SMatthew G. Knepley #if defined(PETSC_USE_64BIT_INDICES) 47235f7792SMatthew G. Knepley inttype = H5T_NATIVE_LLONG; 48235f7792SMatthew G. Knepley #else 49235f7792SMatthew G. Knepley inttype = H5T_NATIVE_INT; 50235f7792SMatthew G. Knepley #endif 51*9566063dSJacob Faibussowitsch PetscCall(PetscObjectGetName((PetscObject)is, &isname)); 52*9566063dSJacob Faibussowitsch PetscCall(PetscViewerHDF5Load(viewer, isname, is->map, inttype, (void**)&ind)); 53*9566063dSJacob Faibussowitsch PetscCall(ISGeneralSetIndices(is, is->map->n, ind, PETSC_OWN_POINTER)); 54235f7792SMatthew G. Knepley PetscFunctionReturn(0); 55235f7792SMatthew G. Knepley } 56235f7792SMatthew G. Knepley #endif 57235f7792SMatthew G. Knepley 58ede126feSBarry Smith PetscErrorCode ISLoad_Binary(IS is, PetscViewer viewer) 59ede126feSBarry Smith { 60f92d9ac3SLisandro Dalcin PetscBool isgeneral,skipHeader; 61f92d9ac3SLisandro Dalcin PetscInt tr[2],rows,N,n,s,*idx; 62f92d9ac3SLisandro Dalcin PetscLayout map; 63ede126feSBarry Smith 64ede126feSBarry Smith PetscFunctionBegin; 65*9566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompare((PetscObject)is,ISGENERAL,&isgeneral)); 66fa7c0e95SVaclav Hapla PetscCheck(isgeneral,PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_INCOMP,"IS must be of type ISGENERAL to load into it"); 67*9566063dSJacob Faibussowitsch PetscCall(PetscViewerSetUp(viewer)); 68*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryGetSkipHeader(viewer,&skipHeader)); 69ede126feSBarry Smith 70*9566063dSJacob Faibussowitsch PetscCall(ISGetLayout(is,&map)); 71*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetSize(map,&N)); 72f92d9ac3SLisandro Dalcin 73f92d9ac3SLisandro Dalcin /* read IS header */ 74f92d9ac3SLisandro Dalcin if (!skipHeader) { 75*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryRead(viewer,tr,2,NULL,PETSC_INT)); 76fa7c0e95SVaclav Hapla PetscCheck(tr[0] == IS_FILE_CLASSID,PetscObjectComm((PetscObject)viewer),PETSC_ERR_FILE_UNEXPECTED,"Not an IS next in file"); 77fa7c0e95SVaclav Hapla PetscCheck(tr[1] >= 0,PetscObjectComm((PetscObject)viewer),PETSC_ERR_FILE_UNEXPECTED,"IS size (%" PetscInt_FMT ") in file is negative",tr[1]); 782c71b3e2SJacob Faibussowitsch PetscCheckFalse(N >= 0 && N != tr[1],PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"IS in file different size (%" PetscInt_FMT ") than input IS (%" PetscInt_FMT ")",tr[1],N); 79f92d9ac3SLisandro Dalcin rows = tr[1]; 80ede126feSBarry Smith } else { 81fa7c0e95SVaclav Hapla PetscCheck(N >= 0,PETSC_COMM_SELF,PETSC_ERR_USER,"IS binary file header was skipped, thus the user must specify the global size of input IS"); 82f92d9ac3SLisandro Dalcin rows = N; 83ede126feSBarry Smith } 84f92d9ac3SLisandro Dalcin 85f92d9ac3SLisandro Dalcin /* set IS size if not already set */ 86*9566063dSJacob Faibussowitsch if (N < 0) PetscCall(PetscLayoutSetSize(map,rows)); 87*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutSetUp(map)); 88f92d9ac3SLisandro Dalcin 89f92d9ac3SLisandro Dalcin /* get IS sizes and check global size */ 90*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetSize(map,&N)); 91*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetLocalSize(map,&n)); 92*9566063dSJacob Faibussowitsch PetscCall(PetscLayoutGetRange(map,&s,NULL)); 93fa7c0e95SVaclav Hapla PetscCheck(N == rows,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"IS in file different size (%" PetscInt_FMT ") than input IS (%" PetscInt_FMT ")",rows,N); 94f92d9ac3SLisandro Dalcin 95f92d9ac3SLisandro Dalcin /* read IS indices */ 96*9566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(n,&idx)); 97*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryReadAll(viewer,idx,n,s,N,PETSC_INT)); 98*9566063dSJacob Faibussowitsch PetscCall(ISGeneralSetIndices(is,n,idx,PETSC_OWN_POINTER)); 99ede126feSBarry Smith PetscFunctionReturn(0); 100ede126feSBarry Smith } 101ede126feSBarry Smith 102235f7792SMatthew G. Knepley PetscErrorCode ISLoad_Default(IS is, PetscViewer viewer) 103235f7792SMatthew G. Knepley { 104f92d9ac3SLisandro Dalcin PetscBool isbinary,ishdf5; 105235f7792SMatthew G. Knepley 106235f7792SMatthew G. Knepley PetscFunctionBegin; 107*9566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary)); 108*9566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5)); 109235f7792SMatthew G. Knepley if (isbinary) { 110*9566063dSJacob Faibussowitsch PetscCall(ISLoad_Binary(is, viewer)); 111235f7792SMatthew G. Knepley } else if (ishdf5) { 112f92d9ac3SLisandro Dalcin #if defined(PETSC_HAVE_HDF5) 113*9566063dSJacob Faibussowitsch PetscCall(ISLoad_HDF5(is, viewer)); 114235f7792SMatthew G. Knepley #endif 115235f7792SMatthew G. Knepley } 116235f7792SMatthew G. Knepley PetscFunctionReturn(0); 117235f7792SMatthew G. Knepley } 118