1 #include <petsc/private/fortranimpl.h> 2 #include <petsc/private/taolinesearchimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define taolinesearchsetobjectiveroutine_ TAOLINESEARCHSETOBJECTIVEROUTINE 6 #define taolinesearchsetgradientroutine_ TAOLINESEARCHSETGRADIENTROUTINE 7 #define taolinesearchsetobjectiveandgradientroutine_ TAOLINESEARCHSETOBJECTIVEANDGRADIENTROUTINE 8 #define taolinesearchsetobjectiveandgtsroutine_ TAOLINESEARCHSETOBJECTIVEANDGTSROUTINE 9 #define taolinesearchview_ TAOLINESEARCHVIEW 10 #define taolinesearchviewfromoptions_ TAOLINESEARCHVIEWFROMOPTIONS 11 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 12 13 #define taolinesearchsetobjectiveroutine_ taolinesearchsetobjectiveroutine 14 #define taolinesearchsetgradientroutine_ taolinesearchsetgradientroutine 15 #define taolinesearchsetobjectiveandgradientroutine_ taolinesearchsetobjectiveandgradientroutine 16 #define taolinesearchsetobjectiveandgtsroutine_ taolinesearchsetobjectiveandgtsroutine 17 #define taolinesearchview_ taolinesearchview 18 #define taolinesearchviewfromoptions_ taolinesearchviewfromoptions 19 #endif 20 21 static int OBJ = 0; 22 static int GRAD = 1; 23 static int OBJGRAD = 2; 24 static int OBJGTS = 3; 25 static size_t NFUNCS = 4; 26 27 static PetscErrorCode ourtaolinesearchobjectiveroutine(TaoLineSearch ls, Vec x, PetscReal *f, void *ctx) 28 { 29 PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, PetscReal *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJ]))(&ls, &x, f, ctx, &ierr)); 30 return PETSC_SUCCESS; 31 } 32 33 static PetscErrorCode ourtaolinesearchgradientroutine(TaoLineSearch ls, Vec x, Vec g, void *ctx) 34 { 35 PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, Vec *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[GRAD]))(&ls, &x, &g, ctx, &ierr)); 36 return PETSC_SUCCESS; 37 } 38 39 static PetscErrorCode ourtaolinesearchobjectiveandgradientroutine(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, void *ctx) 40 { 41 PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJGRAD]))(&ls, &x, f, &g, ctx, &ierr)); 42 return PETSC_SUCCESS; 43 } 44 45 static PetscErrorCode ourtaolinesearchobjectiveandgtsroutine(TaoLineSearch ls, Vec x, Vec s, PetscReal *f, PetscReal *gts, void *ctx) 46 { 47 PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, Vec *, PetscReal *, PetscReal *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJGTS]))(&ls, &x, &s, f, gts, ctx, &ierr)); 48 return PETSC_SUCCESS; 49 } 50 51 PETSC_EXTERN void taolinesearchsetobjectiveroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, PetscReal *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 52 { 53 PetscObjectAllocateFortranPointers(*ls, NFUNCS); 54 if (!func) { 55 *ierr = TaoLineSearchSetObjectiveRoutine(*ls, NULL, ctx); 56 } else { 57 ((PetscObject)*ls)->fortran_func_pointers[OBJ] = (PetscVoidFn *)func; 58 *ierr = TaoLineSearchSetObjectiveRoutine(*ls, ourtaolinesearchobjectiveroutine, ctx); 59 } 60 } 61 62 PETSC_EXTERN void taolinesearchsetgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 63 { 64 PetscObjectAllocateFortranPointers(*ls, NFUNCS); 65 if (!func) { 66 *ierr = TaoLineSearchSetGradientRoutine(*ls, NULL, ctx); 67 } else { 68 ((PetscObject)*ls)->fortran_func_pointers[GRAD] = (PetscVoidFn *)func; 69 *ierr = TaoLineSearchSetGradientRoutine(*ls, ourtaolinesearchgradientroutine, ctx); 70 } 71 } 72 73 PETSC_EXTERN void taolinesearchsetobjectiveandgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 74 { 75 PetscObjectAllocateFortranPointers(*ls, NFUNCS); 76 if (!func) { 77 *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls, NULL, ctx); 78 } else { 79 ((PetscObject)*ls)->fortran_func_pointers[OBJGRAD] = (PetscVoidFn *)func; 80 *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls, ourtaolinesearchobjectiveandgradientroutine, ctx); 81 } 82 } 83 84 PETSC_EXTERN void taolinesearchsetobjectiveandgtsroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, Vec *, PetscReal *, PetscReal *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 85 { 86 PetscObjectAllocateFortranPointers(*ls, NFUNCS); 87 if (!func) { 88 *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls, NULL, ctx); 89 } else { 90 ((PetscObject)*ls)->fortran_func_pointers[OBJGTS] = (PetscVoidFn *)func; 91 *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls, ourtaolinesearchobjectiveandgtsroutine, ctx); 92 } 93 } 94 95 PETSC_EXTERN void taolinesearchview_(TaoLineSearch *ls, PetscViewer *viewer, PetscErrorCode *ierr) 96 { 97 PetscViewer v; 98 PetscPatchDefaultViewers_Fortran(viewer, v); 99 *ierr = TaoLineSearchView(*ls, v); 100 } 101 102 PETSC_EXTERN void taolinesearchviewfromoptions_(TaoLineSearch *ao, PetscObject obj, char *type, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len) 103 { 104 char *t; 105 106 FIXCHAR(type, len, t); 107 CHKFORTRANNULLOBJECT(obj); 108 *ierr = TaoLineSearchViewFromOptions(*ao, obj, t); 109 if (*ierr) return; 110 FREECHAR(type, t); 111 } 112