16dd63270SBarry Smith #include <petsc/private/ftnimpl.h>
2af0996ceSBarry Smith #include <petsc/private/taolinesearchimpl.h>
3a7e14dcfSSatish Balay
4a7e14dcfSSatish Balay #if defined(PETSC_HAVE_FORTRAN_CAPS)
5a7e14dcfSSatish Balay #define taolinesearchsetobjectiveroutine_ TAOLINESEARCHSETOBJECTIVEROUTINE
6a7e14dcfSSatish Balay #define taolinesearchsetgradientroutine_ TAOLINESEARCHSETGRADIENTROUTINE
7a7e14dcfSSatish Balay #define taolinesearchsetobjectiveandgradientroutine_ TAOLINESEARCHSETOBJECTIVEANDGRADIENTROUTINE
8a7e14dcfSSatish Balay #define taolinesearchsetobjectiveandgtsroutine_ TAOLINESEARCHSETOBJECTIVEANDGTSROUTINE
9a7e14dcfSSatish Balay #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
10a7e14dcfSSatish Balay #define taolinesearchsetobjectiveroutine_ taolinesearchsetobjectiveroutine
11a7e14dcfSSatish Balay #define taolinesearchsetgradientroutine_ taolinesearchsetgradientroutine
12a7e14dcfSSatish Balay #define taolinesearchsetobjectiveandgradientroutine_ taolinesearchsetobjectiveandgradientroutine
13a7e14dcfSSatish Balay #define taolinesearchsetobjectiveandgtsroutine_ taolinesearchsetobjectiveandgtsroutine
14a7e14dcfSSatish Balay #endif
15a7e14dcfSSatish Balay
16a7e14dcfSSatish Balay static int OBJ = 0;
17a7e14dcfSSatish Balay static int GRAD = 1;
18a7e14dcfSSatish Balay static int OBJGRAD = 2;
19a7e14dcfSSatish Balay static int OBJGTS = 3;
20e0cd13aeSBarry Smith static size_t NFUNCS = 4;
21a7e14dcfSSatish Balay
ourtaolinesearchobjectiveroutine(TaoLineSearch ls,Vec x,PetscReal * f,PetscCtx ctx)22*2a8381b2SBarry Smith static PetscErrorCode ourtaolinesearchobjectiveroutine(TaoLineSearch ls, Vec x, PetscReal *f, PetscCtx ctx)
23a7e14dcfSSatish Balay {
249566063dSJacob Faibussowitsch PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, PetscReal *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJ]))(&ls, &x, f, ctx, &ierr));
253ba16761SJacob Faibussowitsch return PETSC_SUCCESS;
26a7e14dcfSSatish Balay }
27a7e14dcfSSatish Balay
ourtaolinesearchgradientroutine(TaoLineSearch ls,Vec x,Vec g,PetscCtx ctx)28*2a8381b2SBarry Smith static PetscErrorCode ourtaolinesearchgradientroutine(TaoLineSearch ls, Vec x, Vec g, PetscCtx ctx)
29a7e14dcfSSatish Balay {
309566063dSJacob Faibussowitsch PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, Vec *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[GRAD]))(&ls, &x, &g, ctx, &ierr));
313ba16761SJacob Faibussowitsch return PETSC_SUCCESS;
32a7e14dcfSSatish Balay }
33a7e14dcfSSatish Balay
ourtaolinesearchobjectiveandgradientroutine(TaoLineSearch ls,Vec x,PetscReal * f,Vec g,PetscCtx ctx)34*2a8381b2SBarry Smith static PetscErrorCode ourtaolinesearchobjectiveandgradientroutine(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscCtx ctx)
35a7e14dcfSSatish Balay {
369566063dSJacob Faibussowitsch PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJGRAD]))(&ls, &x, f, &g, ctx, &ierr));
373ba16761SJacob Faibussowitsch return PETSC_SUCCESS;
38a7e14dcfSSatish Balay }
39a7e14dcfSSatish Balay
ourtaolinesearchobjectiveandgtsroutine(TaoLineSearch ls,Vec x,Vec s,PetscReal * f,PetscReal * gts,PetscCtx ctx)40*2a8381b2SBarry Smith static PetscErrorCode ourtaolinesearchobjectiveandgtsroutine(TaoLineSearch ls, Vec x, Vec s, PetscReal *f, PetscReal *gts, PetscCtx ctx)
41a7e14dcfSSatish Balay {
429566063dSJacob Faibussowitsch PetscCallFortranVoidFunction((*(void (*)(TaoLineSearch *, Vec *, Vec *, PetscReal *, PetscReal *, void *, PetscErrorCode *))(((PetscObject)ls)->fortran_func_pointers[OBJGTS]))(&ls, &x, &s, f, gts, ctx, &ierr));
433ba16761SJacob Faibussowitsch return PETSC_SUCCESS;
44a7e14dcfSSatish Balay }
45a7e14dcfSSatish Balay
taolinesearchsetobjectiveroutine_(TaoLineSearch * ls,void (* func)(TaoLineSearch *,Vec *,PetscReal *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)46*2a8381b2SBarry Smith PETSC_EXTERN void taolinesearchsetobjectiveroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, PetscReal *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
47a7e14dcfSSatish Balay {
48a7e14dcfSSatish Balay PetscObjectAllocateFortranPointers(*ls, NFUNCS);
49a7e14dcfSSatish Balay if (!func) {
50dfef5ea7SSatish Balay *ierr = TaoLineSearchSetObjectiveRoutine(*ls, NULL, ctx);
51a7e14dcfSSatish Balay } else {
525ebfa9e9SBarry Smith ((PetscObject)*ls)->fortran_func_pointers[OBJ] = (PetscFortranCallbackFn *)func;
53a7e14dcfSSatish Balay *ierr = TaoLineSearchSetObjectiveRoutine(*ls, ourtaolinesearchobjectiveroutine, ctx);
54a7e14dcfSSatish Balay }
55a7e14dcfSSatish Balay }
56a7e14dcfSSatish Balay
taolinesearchsetgradientroutine_(TaoLineSearch * ls,void (* func)(TaoLineSearch *,Vec *,Vec *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)57*2a8381b2SBarry Smith PETSC_EXTERN void taolinesearchsetgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, Vec *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
58a7e14dcfSSatish Balay {
59a7e14dcfSSatish Balay PetscObjectAllocateFortranPointers(*ls, NFUNCS);
60a7e14dcfSSatish Balay if (!func) {
61dfef5ea7SSatish Balay *ierr = TaoLineSearchSetGradientRoutine(*ls, NULL, ctx);
62a7e14dcfSSatish Balay } else {
635ebfa9e9SBarry Smith ((PetscObject)*ls)->fortran_func_pointers[GRAD] = (PetscFortranCallbackFn *)func;
64a7e14dcfSSatish Balay *ierr = TaoLineSearchSetGradientRoutine(*ls, ourtaolinesearchgradientroutine, ctx);
65a7e14dcfSSatish Balay }
66a7e14dcfSSatish Balay }
67a7e14dcfSSatish Balay
taolinesearchsetobjectiveandgradientroutine_(TaoLineSearch * ls,void (* func)(TaoLineSearch *,Vec *,PetscReal *,Vec *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)68*2a8381b2SBarry Smith PETSC_EXTERN void taolinesearchsetobjectiveandgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
69a7e14dcfSSatish Balay {
70a7e14dcfSSatish Balay PetscObjectAllocateFortranPointers(*ls, NFUNCS);
71a7e14dcfSSatish Balay if (!func) {
72dfef5ea7SSatish Balay *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls, NULL, ctx);
73a7e14dcfSSatish Balay } else {
745ebfa9e9SBarry Smith ((PetscObject)*ls)->fortran_func_pointers[OBJGRAD] = (PetscFortranCallbackFn *)func;
75a7e14dcfSSatish Balay *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls, ourtaolinesearchobjectiveandgradientroutine, ctx);
76a7e14dcfSSatish Balay }
77a7e14dcfSSatish Balay }
78a7e14dcfSSatish Balay
taolinesearchsetobjectiveandgtsroutine_(TaoLineSearch * ls,void (* func)(TaoLineSearch *,Vec *,Vec *,PetscReal *,PetscReal *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)79*2a8381b2SBarry Smith PETSC_EXTERN void taolinesearchsetobjectiveandgtsroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch *, Vec *, Vec *, PetscReal *, PetscReal *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
80a7e14dcfSSatish Balay {
81a7e14dcfSSatish Balay PetscObjectAllocateFortranPointers(*ls, NFUNCS);
82a7e14dcfSSatish Balay if (!func) {
83dfef5ea7SSatish Balay *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls, NULL, ctx);
84a7e14dcfSSatish Balay } else {
855ebfa9e9SBarry Smith ((PetscObject)*ls)->fortran_func_pointers[OBJGTS] = (PetscFortranCallbackFn *)func;
86a7e14dcfSSatish Balay *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls, ourtaolinesearchobjectiveandgtsroutine, ctx);
87a7e14dcfSSatish Balay }
88a7e14dcfSSatish Balay }
89