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 19f6291634SJed Brown static PetscErrorCode PetscFortranCallbackFinalize(void) 20f6291634SJed Brown { 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; 27*9566063dSJacob Faibussowitsch PetscCall(PetscFree(link->type_name)); 28*9566063dSJacob Faibussowitsch PetscCall(PetscFree(link)); 29f6291634SJed Brown } 30f6291634SJed Brown } 31*9566063dSJacob Faibussowitsch PetscCall(PetscFree(_classbase)); 32f6291634SJed Brown _maxclassid = PETSC_SMALLEST_CLASSID; 33f6291634SJed Brown PetscFunctionReturn(0); 34f6291634SJed Brown } 35f6291634SJed Brown 36de6d466bSJed Brown /*@C 37f6291634SJed Brown PetscFortranCallbackRegister - register a type+subtype callback 38f6291634SJed Brown 39f6291634SJed Brown Not Collective 40f6291634SJed Brown 414165533cSJose E. Roman Input Parameters: 42f6291634SJed Brown + classid - ID of class on which to register callback 430298fd71SBarry Smith - subtype - subtype string, or NULL for class ids 44f6291634SJed Brown 454165533cSJose E. Roman Output Parameter: 46f6291634SJed Brown . id - callback id 47f6291634SJed Brown 48f6291634SJed Brown Level: developer 49f6291634SJed Brown 50f6291634SJed Brown .seealso: PetscFortranCallbackGetSizes() 51f6291634SJed Brown @*/ 52f6291634SJed Brown PetscErrorCode PetscFortranCallbackRegister(PetscClassId classid,const char *subtype,PetscFortranCallbackId *id) 53f6291634SJed Brown { 54f6291634SJed Brown FortranCallbackBase *base; 55f6291634SJed Brown FortranCallbackLink link; 56f6291634SJed Brown 57f6291634SJed Brown PetscFunctionBegin; 583ca90d2dSJacob Faibussowitsch if (subtype) PetscValidCharPointer(subtype,2); 593ca90d2dSJacob Faibussowitsch PetscValidPointer(id,3); 602c71b3e2SJacob Faibussowitsch PetscCheckFalse(classid < PETSC_SMALLEST_CLASSID || PETSC_LARGEST_CLASSID < classid,PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"ClassId %d corrupt",classid); 61f6291634SJed Brown *id = 0; 62f6291634SJed Brown if (classid >= _maxclassid) { 63f6291634SJed Brown PetscClassId newmax = PETSC_SMALLEST_CLASSID + 2*(PETSC_LARGEST_CLASSID-PETSC_SMALLEST_CLASSID); 64f6291634SJed Brown FortranCallbackBase *newbase; 65*9566063dSJacob Faibussowitsch if (!_classbase) PetscCall(PetscRegisterFinalize(PetscFortranCallbackFinalize)); 66*9566063dSJacob Faibussowitsch PetscCall(PetscCalloc1(newmax-PETSC_SMALLEST_CLASSID,&newbase)); 67*9566063dSJacob Faibussowitsch PetscCall(PetscArraycpy(newbase,_classbase,_maxclassid-PETSC_SMALLEST_CLASSID)); 68*9566063dSJacob Faibussowitsch PetscCall(PetscFree(_classbase)); 69a297a907SKarl Rupp 70f6291634SJed Brown _classbase = newbase; 71f6291634SJed Brown _maxclassid = newmax; 72f6291634SJed Brown } 73f6291634SJed Brown base = &_classbase[classid-PETSC_SMALLEST_CLASSID]; 74a297a907SKarl Rupp if (!subtype) *id = PETSC_SMALLEST_FORTRAN_CALLBACK + base->basecount++; 75a297a907SKarl Rupp else { 76f6291634SJed Brown for (link=base->subtypes; link; link=link->next) { /* look for either both NULL or matching values (implies both non-NULL) */ 77f6291634SJed Brown PetscBool match; 78*9566063dSJacob Faibussowitsch PetscCall(PetscStrcmp(subtype,link->type_name,&match)); 79f6291634SJed Brown if (match) { /* base type or matching subtype */ 80f6291634SJed Brown goto found; 81f6291634SJed Brown } 82f6291634SJed Brown } 83f6291634SJed Brown /* Not found. Create node and prepend to class' subtype list */ 84*9566063dSJacob Faibussowitsch PetscCall(PetscNew(&link)); 85*9566063dSJacob Faibussowitsch PetscCall(PetscStrallocpy(subtype,&link->type_name)); 86a297a907SKarl Rupp 87f6291634SJed Brown link->max = PETSC_SMALLEST_FORTRAN_CALLBACK; 88f6291634SJed Brown link->next = base->subtypes; 89f6291634SJed Brown base->subtypes = link; 90a297a907SKarl Rupp 91f6291634SJed Brown found: 92f6291634SJed Brown *id = link->max++; 93a297a907SKarl Rupp 94f6291634SJed Brown base->maxsubtypecount = PetscMax(base->maxsubtypecount,link->max-PETSC_SMALLEST_FORTRAN_CALLBACK); 95f6291634SJed Brown } 96f6291634SJed Brown PetscFunctionReturn(0); 97f6291634SJed Brown } 98f6291634SJed Brown 99de6d466bSJed Brown /*@C 100f6291634SJed Brown PetscFortranCallbackGetSizes - get sizes of class and subtype pointer arrays 101f6291634SJed Brown 102f6291634SJed Brown Collective 103f6291634SJed Brown 1044165533cSJose E. Roman Input Parameter: 105f6291634SJed Brown . classid - class Id 106f6291634SJed Brown 1074165533cSJose E. Roman Output Parameters: 108f6291634SJed Brown + numbase - number of registered class callbacks 109f6291634SJed Brown - numsubtype - max number of registered subtype callbacks 110f6291634SJed Brown 111f6291634SJed Brown Level: developer 112f6291634SJed Brown 113f6291634SJed Brown .seealso: PetscFortranCallbackRegister() 114f6291634SJed Brown @*/ 115e0cd13aeSBarry Smith PetscErrorCode PetscFortranCallbackGetSizes(PetscClassId classid,PetscFortranCallbackId *numbase,PetscFortranCallbackId *numsubtype) 116f6291634SJed Brown { 117f6291634SJed Brown PetscFunctionBegin; 1185f80ce2aSJacob Faibussowitsch PetscValidPointer(numbase,2); 1195f80ce2aSJacob Faibussowitsch PetscValidPointer(numsubtype,3); 120f6291634SJed Brown if (classid < _maxclassid) { 121f6291634SJed Brown FortranCallbackBase *base = &_classbase[classid-PETSC_SMALLEST_CLASSID]; 122f6291634SJed Brown *numbase = base->basecount; 123f6291634SJed Brown *numsubtype = base->maxsubtypecount; 124f6291634SJed Brown } else { /* nothing registered */ 125f6291634SJed Brown *numbase = 0; 126f6291634SJed Brown *numsubtype = 0; 127f6291634SJed Brown } 128f6291634SJed Brown PetscFunctionReturn(0); 129f6291634SJed Brown } 130