On Tue, Nov 28, 2023 at 5:23 PM Brandon Denton <blden...@buffalo.edu> wrote:
> Good Afternoon, > > Ok. I think I've developed a solution that meets everyone's interest. > Please review the below and let me know if you have any questions or > comments. > > Quick background. The objective of this function is to provide user > friendly access to geometric data associated with a DMPlex with a geometry > CAD model attached. This function will replace the originally conceived > PetscObjectQuery() method. > Yes, that is a good start. We might eventually want an object, but this is great to experiment with. I like it. Thanks, Matt > *Proposed Function Signature:* > DMPlexGetGeomCntrlPntAndWeightData(DM dm, PetscHMapI *cpHashTable, > PetscInt *cpCoordDataLength, PetscScalar **cpCoordData, PetscInt > *maxNumEquiv, Mat *cpEquiv, PetscHMapI *wHashTable, PetscInt *wDataLength, > PetscScalar **wData) > > The only input it the DM. All others are outputs. > > *Proposed Call* > PetscScalar *cpCoordData, *wData; > PetscInt cpCoordDataLength, wDataLength, maxNumEquiv; > PetscHMapI cpHashTable, wHashTable; > Mat cpEquiv; > > PetscCall(DMPlexGetGeomCntrlPntAndWeightData(dmNozzle, &cpHashTable, > &cpCoordDataLength, &cpCoordData, &maxNumEquiv, &cpEquiv, &wHashTable, > &wDataLength, &wData)); > > *Proposed Function Added to plexegads.c* > DMPlexGetGeomCntrlPntAndWeightData(DM dm, PetscHMapI *cpHashTable, > PetscInt *cpCoordDataLength, PetscScalar **cpCoordData, PetscInt > *maxNumEquiv, Mat *cpEquiv, PetscHMapI *wHashTable, PetscInt *wDataLength, > PetscScalar **wData) > { > PetscFunctionBeginHot; > #ifdef PETSC_HAVE_EGADS > PetscContainer modelObj, cpHashTableObj, wHashTableObj, cpCoordDataObj, > wDataObj, cpCoordDataLengthObj, wDataLengthObj, cpEquivObj, maxNumRelateObj; > PetscInt *cpCoordDataLengthPtr, *wDataLengthPtr, *maxNumEquivPtr; > PetscHMapI cpHashTableTemp, wHashTableTemp; > PetscScalar *cpCoordDataTemp, *wDataTemp; > Mat cpEquivTemp; > > /* Determine which type of EGADS model is attached to the DM */ > PetscCall(PetscObjectQuery((PetscObject) dm, "EGADS Model", (PetscObject > *) &modelObj)); > if (!modelObj) { > PetscCall(PetscObjectQuery((PetscObject) dm, "EGADSlite Model", > (PetscObject *) &modelObj)); > } > > if (!modelObj) {PetscFunctionReturn(PETSC_SUCCESS);} > > // Look to see if DM has Container for Geometry Control Point Data > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Hash Table", > (PetscObject *)&cpHashTableObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Coordinates", > (PetscObject *)&cpCoordDataObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Coordinate > Data Length", (PetscObject *)&cpCoordDataLengthObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Weights Hash > Table", (PetscObject *)&wHashTableObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Weight Data", > (PetscObject *)&wDataObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Weight Data > Length", (PetscObject *)&wDataLengthObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Control Point Equivalancy > Matrix", (PetscObject *)&cpEquivObj)); > PetscCall(PetscObjectQuery((PetscObject)dm, "Maximum Number Control > Point Equivalency", (PetscObject *)&maxNumRelateObj)); > > // Get attached EGADS model Control Point and Weights Hash Tables and > Data Arrays (pointer) > PetscCall(PetscContainerGetPointer(cpHashTableObj, (void > **)&cpHashTableTemp)); > PetscCall(PetscContainerGetPointer(cpCoordDataObj, (void > **)&cpCoordDataTemp)); > PetscCall(PetscContainerGetPointer(cpCoordDataLengthObj, (void > **)&cpCoordDataLengthPtr)); > PetscCall(PetscContainerGetPointer(wHashTableObj, (void > **)&wHashTableTemp)); > PetscCall(PetscContainerGetPointer(wDataObj, (void **)&wDataTemp)); > PetscCall(PetscContainerGetPointer(wDataLengthObj, (void > **)&wDataLengthPtr)); > PetscCall(PetscContainerGetPointer(cpEquivObj, (void **)&cpEquivTemp)); > PetscCall(PetscContainerGetPointer(maxNumRelateObj, (void > **)&maxNumEquivPtr)); > > *cpCoordDataLength = *cpCoordDataLengthPtr; > *wDataLength = *wDataLengthPtr; > *maxNumEquiv = *maxNumEquivPtr; > *cpEquiv = cpEquivTemp; > *cpHashTable = cpHashTableTemp; > *wHashTable = wHashTableTemp; > *cpCoordData = cpCoordDataTemp; > *wData = wDataTemp; > #endif > > PetscFunctionReturn(PETSC_SUCCESS); > } > > ------------------------------ > *From:* Matthew Knepley <knep...@gmail.com> > *Sent:* Tuesday, November 28, 2023 1:26 PM > *To:* Barry Smith <bsm...@petsc.dev> > *Cc:* Brandon Denton <blden...@buffalo.edu>; petsc-dev@mcs.anl.gov < > petsc-dev@mcs.anl.gov> > *Subject:* Re: [petsc-dev] PetscHMapI passing > > On Mon, Nov 27, 2023 at 1:00 PM Barry Smith <bsm...@petsc.dev> wrote: > > > What errors do you get? > > I am trying to develop a function where the user can specify an empty > PetscHMapI variable and the function will retrieve the information attached > to the DMPlex. However, when I try something like ... > > PETSC_EXTERN PetscErrorCode DMPlexGetGeomInfo(DM dm, PetscHMapI hashTable) > > > I don't understand the design of specifying and empty one. Why not have > > PETSC_EXTERN PetscErrorCode DMPlexGetGeomInfo(DM dm, PetscHMapI > *hashTable) > > and return the one you have inside your container? > > > This sounds like the right design. Does that work? > > Thanks, > > Matt > > > On Nov 27, 2023, at 9:11 AM, Brandon Denton via petsc-dev < > petsc-dev@mcs.anl.gov> wrote: > > Good Morning, > > I am trying to develop a function in PETSc that allows for the > retrieval of a PetscHMapI contained in a PetscContainer attached to a > DMPlex. The creation of this Hash Table and its attachment to the DMPlex is > automated as part of the CAD/Discretization integration I've been > developing. Currently, the only way to get this Hash Table is through the > PetscObjectQuery() functionality. In an effort to make it easier for users > to access and use this information, I am trying to develop a function where > the user can specify an empty PetscHMapI variable and the function will > retrieve the information attached to the DMPlex. However, when I try > something like ... > > PETSC_EXTERN PetscErrorCode DMPlexGetGeomInfo(DM dm, PetscHMapI hashTable) > > I get a ton of compilation errors associated with the PetscHMapI > designation. Is there a way around this? > > Thank you. > Brandon > > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/~knepley/> > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>