1 #include <petsc/private/ftnimpl.h> 2 #include <petsc/private/taoimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define taoadmmsetmisfitobjectiveandgradientroutine_ TAOADMMSETMISFITOBJECTIVEANDGRADIENTROUTINE 6 #define taoadmmsetmisfithessianroutine_ TAOADMMSETMISFITHESSIANROUTINE 7 #define taoadmmsetmisfitconstraintjacobian_ TAOADMMSETMISFITCONSTRAINTJACOBIAN 8 #define taoadmmsetregularizerobjectiveandgradientroutine_ TAOADMMSETREGULARIZEROBJECTIVEANDGRADIENTROUTINE 9 #define taoadmmsetregularizerhessianroutine_ TAOADMMSETREGULARIZERHESSIANROUTINE 10 #define taoadmmsetregularizerconstraintjacobian_ TAOADMMSETREGULARIZERCONSTRAINTJACOBIAN 11 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 12 #define taoadmmsetmisfitobjectiveandgradientroutine_ taoadmmsetmisfitobjectiveandgradientroutine 13 #define taoadmmsetmisfithessianroutine_ taoadmmsetmisfithessianroutine 14 #define taoadmmsetmisfitconstraintjacobian_ taoadmmsetmisfitconstraintjacobian 15 #define taoadmmsetregularizerobjectiveandgradientroutine_ taoadmmsetregularizerobjectiveandgradientroutine 16 #define taoadmmsetregularizerhessianroutine_ taoadmmsetregularizerhessianroutine 17 #define taoadmmsetregularizerconstraintjacobian_ taoadmmsetregularizerconstraintjacobian 18 #endif 19 20 static struct { 21 PetscFortranCallbackId misfitobjgrad; 22 PetscFortranCallbackId misfithess; 23 PetscFortranCallbackId misfitjacobian; 24 PetscFortranCallbackId regobjgrad; 25 PetscFortranCallbackId reghess; 26 PetscFortranCallbackId regjacobian; 27 } _cb; 28 29 static PetscErrorCode ourtaoadmmmisfitobjgradroutine(Tao tao, Vec x, PetscReal *f, Vec g, void *ctx) 30 { 31 PetscObjectUseFortranCallback(tao, _cb.misfitobjgrad, (Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), (&tao, &x, f, &g, _ctx, &ierr)); 32 } 33 34 static PetscErrorCode ourtaoadmmmisfithessroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx) 35 { 36 PetscObjectUseFortranCallback(tao, _cb.misfithess, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &H, &Hpre, _ctx, &ierr)); 37 } 38 39 static PetscErrorCode ourtaoadmmmisfitconstraintjacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx) 40 { 41 PetscObjectUseFortranCallback(tao, _cb.misfitjacobian, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &J, &Jpre, _ctx, &ierr)); 42 } 43 44 static PetscErrorCode ourtaoadmmregularizerobjgradroutine(Tao tao, Vec x, PetscReal *f, Vec g, void *ctx) 45 { 46 PetscObjectUseFortranCallback(tao, _cb.regobjgrad, (Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), (&tao, &x, f, &g, _ctx, &ierr)); 47 } 48 49 static PetscErrorCode ourtaoadmmregularizerhessroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx) 50 { 51 PetscObjectUseFortranCallback(tao, _cb.reghess, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &H, &Hpre, _ctx, &ierr)); 52 } 53 54 static PetscErrorCode ourtaoadmmregularizerconstraintjacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx) 55 { 56 PetscObjectUseFortranCallback(tao, _cb.regjacobian, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &J, &Jpre, _ctx, &ierr)); 57 } 58 59 PETSC_EXTERN void taoadmmsetmisfitobjectiveandgradientroutine_(Tao *tao, void (*func)(Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 60 { 61 CHKFORTRANNULLFUNCTION(func); 62 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitobjgrad, (PetscFortranCallbackFn *)func, ctx); 63 if (!*ierr) *ierr = TaoADMMSetMisfitObjectiveAndGradientRoutine(*tao, ourtaoadmmmisfitobjgradroutine, ctx); 64 } 65 66 PETSC_EXTERN void taoadmmsetmisfithessianroutine_(Tao *tao, Mat *H, Mat *Hpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 67 { 68 CHKFORTRANNULLFUNCTION(func); 69 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfithess, (PetscFortranCallbackFn *)func, ctx); 70 if (!*ierr) *ierr = TaoADMMSetMisfitHessianRoutine(*tao, *H, *Hpre, ourtaoadmmmisfithessroutine, ctx); 71 } 72 73 PETSC_EXTERN void taoadmmsetmisfitconstraintjacobian_(Tao *tao, Mat *J, Mat *Jpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 74 { 75 CHKFORTRANNULLFUNCTION(func); 76 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitjacobian, (PetscFortranCallbackFn *)func, ctx); 77 if (!*ierr) *ierr = TaoADMMSetMisfitConstraintJacobian(*tao, *J, *Jpre, ourtaoadmmmisfitconstraintjacobian, ctx); 78 } 79 80 PETSC_EXTERN void taoadmmsetregularizerobjectiveandgradientroutine_(Tao *tao, void (*func)(Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 81 { 82 CHKFORTRANNULLFUNCTION(func); 83 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.regobjgrad, (PetscFortranCallbackFn *)func, ctx); 84 if (!*ierr) *ierr = TaoADMMSetRegularizerObjectiveAndGradientRoutine(*tao, ourtaoadmmregularizerobjgradroutine, ctx); 85 } 86 87 PETSC_EXTERN void taoadmmsetregularizerhessianroutine_(Tao *tao, Mat *H, Mat *Hpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 88 { 89 CHKFORTRANNULLFUNCTION(func); 90 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.reghess, (PetscFortranCallbackFn *)func, ctx); 91 if (!*ierr) *ierr = TaoADMMSetRegularizerHessianRoutine(*tao, *H, *Hpre, ourtaoadmmregularizerhessroutine, ctx); 92 } 93 94 PETSC_EXTERN void taoadmmsetregularizerconstraintjacobian_(Tao *tao, Mat *J, Mat *Jpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 95 { 96 CHKFORTRANNULLFUNCTION(func); 97 *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitjacobian, (PetscFortranCallbackFn *)func, ctx); 98 if (!*ierr) *ierr = TaoADMMSetRegularizerConstraintJacobian(*tao, *J, *Jpre, ourtaoadmmregularizerconstraintjacobian, ctx); 99 } 100