1af0996ceSBarry Smith #include <petsc/private/petscimpl.h> /*I "petscsys.h" I*/ 2f6291634SJed Brown 3f6291634SJed Brown typedef struct _FortranCallbackLink *FortranCallbackLink; 4f6291634SJed Brown struct _FortranCallbackLink { 5f6291634SJed Brown char *type_name; 6f6291634SJed Brown PetscFortranCallbackId max; 7f6291634SJed Brown FortranCallbackLink next; 8f6291634SJed Brown }; 9f6291634SJed Brown 10f6291634SJed Brown typedef struct { 11e0cd13aeSBarry Smith PetscFortranCallbackId basecount; 12e0cd13aeSBarry Smith PetscFortranCallbackId maxsubtypecount; 13f6291634SJed Brown FortranCallbackLink subtypes; 14f6291634SJed Brown } FortranCallbackBase; 15f6291634SJed Brown 16f6291634SJed Brown static FortranCallbackBase *_classbase; 17f6291634SJed Brown static PetscClassId _maxclassid = PETSC_SMALLEST_CLASSID; 18f6291634SJed Brown 19*d71ae5a4SJacob Faibussowitsch static PetscErrorCode PetscFortranCallbackFinalize(void) 20*d71ae5a4SJacob Faibussowitsch { 21f6291634SJed Brown PetscFunctionBegin; 225f80ce2aSJacob Faibussowitsch for (PetscInt i = PETSC_SMALLEST_CLASSID; i < _maxclassid; i++) { 23f6291634SJed Brown FortranCallbackBase *base = &_classbase[i - PETSC_SMALLEST_CLASSID]; 24f6291634SJed Brown FortranCallbackLink next, link = base->subtypes; 25f6291634SJed Brown for (; link; link = next) { 26f6291634SJed Brown next = link->next; 279566063dSJacob Faibussowitsch PetscCall(PetscFree(link->type_name)); 289566063dSJacob Faibussowitsch PetscCall(PetscFree(link)); 29f6291634SJed Brown } 30f6291634SJed Brown } 319566063dSJacob Faibussowitsch PetscCall(PetscFree(_classbase)); 32f6291634SJed Brown _maxclassid = PETSC_SMALLEST_CLASSID; 33f6291634SJed Brown PetscFunctionReturn(0); 34f6291634SJed Brown } 35f6291634SJed Brown 36de6d466bSJed Brown /*@C 37811af0c4SBarry Smith PetscFortranCallbackRegister - register a type+subtype callback. This is used by the PETSc Fortran interface to allow the use of user Fortran functions as arguments 38811af0c4SBarry Smith to PETSc functions that take function pointers 39f6291634SJed Brown 40f6291634SJed Brown Not Collective 41f6291634SJed Brown 424165533cSJose E. Roman Input Parameters: 43f6291634SJed Brown + classid - ID of class on which to register callback 440298fd71SBarry Smith - subtype - subtype string, or NULL for class ids 45f6291634SJed Brown 464165533cSJose E. Roman Output Parameter: 47f6291634SJed Brown . id - callback id 48f6291634SJed Brown 49f6291634SJed Brown Level: developer 50f6291634SJed Brown 51db781477SPatrick Sanan .seealso: `PetscFortranCallbackGetSizes()` 52f6291634SJed Brown @*/ 53*d71ae5a4SJacob Faibussowitsch PetscErrorCode PetscFortranCallbackRegister(PetscClassId classid, const char *subtype, PetscFortranCallbackId *id) 54*d71ae5a4SJacob Faibussowitsch { 55f6291634SJed Brown FortranCallbackBase *base; 56f6291634SJed Brown FortranCallbackLink link; 57f6291634SJed Brown 58f6291634SJed Brown PetscFunctionBegin; 593ca90d2dSJacob Faibussowitsch if (subtype) PetscValidCharPointer(subtype, 2); 603ca90d2dSJacob Faibussowitsch PetscValidPointer(id, 3); 61cc73adaaSBarry Smith PetscCheck(classid >= PETSC_SMALLEST_CLASSID && classid <= PETSC_LARGEST_CLASSID, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "ClassId %d corrupt", classid); 62f6291634SJed Brown *id = 0; 63f6291634SJed Brown if (classid >= _maxclassid) { 64f6291634SJed Brown PetscClassId newmax = PETSC_SMALLEST_CLASSID + 2 * (PETSC_LARGEST_CLASSID - PETSC_SMALLEST_CLASSID); 65f6291634SJed Brown FortranCallbackBase *newbase; 669566063dSJacob Faibussowitsch if (!_classbase) PetscCall(PetscRegisterFinalize(PetscFortranCallbackFinalize)); 679566063dSJacob Faibussowitsch PetscCall(PetscCalloc1(newmax - PETSC_SMALLEST_CLASSID, &newbase)); 689566063dSJacob Faibussowitsch PetscCall(PetscArraycpy(newbase, _classbase, _maxclassid - PETSC_SMALLEST_CLASSID)); 699566063dSJacob Faibussowitsch PetscCall(PetscFree(_classbase)); 70a297a907SKarl Rupp 71f6291634SJed Brown _classbase = newbase; 72f6291634SJed Brown _maxclassid = newmax; 73f6291634SJed Brown } 74f6291634SJed Brown base = &_classbase[classid - PETSC_SMALLEST_CLASSID]; 75a297a907SKarl Rupp if (!subtype) *id = PETSC_SMALLEST_FORTRAN_CALLBACK + base->basecount++; 76a297a907SKarl Rupp else { 77f6291634SJed Brown for (link = base->subtypes; link; link = link->next) { /* look for either both NULL or matching values (implies both non-NULL) */ 78f6291634SJed Brown PetscBool match; 799566063dSJacob Faibussowitsch PetscCall(PetscStrcmp(subtype, link->type_name, &match)); 80f6291634SJed Brown if (match) { /* base type or matching subtype */ 81f6291634SJed Brown goto found; 82f6291634SJed Brown } 83f6291634SJed Brown } 84f6291634SJed Brown /* Not found. Create node and prepend to class' subtype list */ 859566063dSJacob Faibussowitsch PetscCall(PetscNew(&link)); 869566063dSJacob Faibussowitsch PetscCall(PetscStrallocpy(subtype, &link->type_name)); 87a297a907SKarl Rupp 88f6291634SJed Brown link->max = PETSC_SMALLEST_FORTRAN_CALLBACK; 89f6291634SJed Brown link->next = base->subtypes; 90f6291634SJed Brown base->subtypes = link; 91a297a907SKarl Rupp 92f6291634SJed Brown found: 93f6291634SJed Brown *id = link->max++; 94a297a907SKarl Rupp 95f6291634SJed Brown base->maxsubtypecount = PetscMax(base->maxsubtypecount, link->max - PETSC_SMALLEST_FORTRAN_CALLBACK); 96f6291634SJed Brown } 97f6291634SJed Brown PetscFunctionReturn(0); 98f6291634SJed Brown } 99f6291634SJed Brown 100de6d466bSJed Brown /*@C 101f6291634SJed Brown PetscFortranCallbackGetSizes - get sizes of class and subtype pointer arrays 102f6291634SJed Brown 103f6291634SJed Brown Collective 104f6291634SJed Brown 1054165533cSJose E. Roman Input Parameter: 106f6291634SJed Brown . classid - class Id 107f6291634SJed Brown 1084165533cSJose E. Roman Output Parameters: 109f6291634SJed Brown + numbase - number of registered class callbacks 110f6291634SJed Brown - numsubtype - max number of registered subtype callbacks 111f6291634SJed Brown 112f6291634SJed Brown Level: developer 113f6291634SJed Brown 114db781477SPatrick Sanan .seealso: `PetscFortranCallbackRegister()` 115f6291634SJed Brown @*/ 116*d71ae5a4SJacob Faibussowitsch PetscErrorCode PetscFortranCallbackGetSizes(PetscClassId classid, PetscFortranCallbackId *numbase, PetscFortranCallbackId *numsubtype) 117*d71ae5a4SJacob Faibussowitsch { 118f6291634SJed Brown PetscFunctionBegin; 1195f80ce2aSJacob Faibussowitsch PetscValidPointer(numbase, 2); 1205f80ce2aSJacob Faibussowitsch PetscValidPointer(numsubtype, 3); 121f6291634SJed Brown if (classid < _maxclassid) { 122f6291634SJed Brown FortranCallbackBase *base = &_classbase[classid - PETSC_SMALLEST_CLASSID]; 123f6291634SJed Brown *numbase = base->basecount; 124f6291634SJed Brown *numsubtype = base->maxsubtypecount; 125f6291634SJed Brown } else { /* nothing registered */ 126f6291634SJed Brown *numbase = 0; 127f6291634SJed Brown *numsubtype = 0; 128f6291634SJed Brown } 129f6291634SJed Brown PetscFunctionReturn(0); 130f6291634SJed Brown } 131