Actual source code: zadmmf.c

  1: #include <petsc/private/fortranimpl.h>
  2: #include <petsc/private/f90impl.h>
  3: #include <petsc/private/taoimpl.h>

  5: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  6:   #define taoadmmsetmisfitobjectiveandgradientroutine_      TAOADMMSETMISFITOBJECTIVEANDGRADIENTROUTINE
  7:   #define taoadmmsetmisfithessianroutine_                   TAOADMMSETMISFITHESSIANROUTINE
  8:   #define taoadmmsetmisfitconstraintjacobian_               TAOADMMSETMISFITCONSTRAINTJACOBIAN
  9:   #define taoadmmsetregularizerobjectiveandgradientroutine_ TAOADMMSETREGULARIZEROBJECTIVEANDGRADIENTROUTINE
 10:   #define taoadmmsetregularizerhessianroutine_              TAOADMMSETREGULARIZERHESSIANROUTINE
 11:   #define taoadmmsetregularizerconstraintjacobian_          TAOADMMSETREGULARIZERCONSTRAINTJACOBIAN
 12: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 13:   #define taoadmmsetmisfitobjectiveandgradientroutine_      taoadmmsetmisfitobjectiveandgradientroutine
 14:   #define taoadmmsetmisfithessianroutine_                   taoadmmsetmisfithessianroutine
 15:   #define taoadmmsetmisfitconstraintjacobian_               taoadmmsetmisfitconstraintjacobian
 16:   #define taoadmmsetregularizerobjectiveandgradientroutine_ taoadmmsetregularizerobjectiveandgradientroutine
 17:   #define taoadmmsetregularizerhessianroutine_              taoadmmsetregularizerhessianroutine
 18:   #define taoadmmsetregularizerconstraintjacobian_          taoadmmsetregularizerconstraintjacobian
 19: #endif

 21: static struct {
 22:   PetscFortranCallbackId misfitobjgrad;
 23:   PetscFortranCallbackId misfithess;
 24:   PetscFortranCallbackId misfitjacobian;
 25:   PetscFortranCallbackId regobjgrad;
 26:   PetscFortranCallbackId reghess;
 27:   PetscFortranCallbackId regjacobian;
 28: } _cb;

 30: static PetscErrorCode ourtaoadmmmisfitobjgradroutine(Tao tao, Vec x, PetscReal *f, Vec g, void *ctx)
 31: {
 32:   PetscObjectUseFortranCallback(tao, _cb.misfitobjgrad, (Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), (&tao, &x, f, &g, _ctx, &ierr));
 33: }

 35: static PetscErrorCode ourtaoadmmmisfithessroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx)
 36: {
 37:   PetscObjectUseFortranCallback(tao, _cb.misfithess, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &H, &Hpre, _ctx, &ierr));
 38: }

 40: static PetscErrorCode ourtaoadmmmisfitconstraintjacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
 41: {
 42:   PetscObjectUseFortranCallback(tao, _cb.misfitjacobian, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &J, &Jpre, _ctx, &ierr));
 43: }

 45: static PetscErrorCode ourtaoadmmregularizerobjgradroutine(Tao tao, Vec x, PetscReal *f, Vec g, void *ctx)
 46: {
 47:   PetscObjectUseFortranCallback(tao, _cb.regobjgrad, (Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), (&tao, &x, f, &g, _ctx, &ierr));
 48: }

 50: static PetscErrorCode ourtaoadmmregularizerhessroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx)
 51: {
 52:   PetscObjectUseFortranCallback(tao, _cb.reghess, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &H, &Hpre, _ctx, &ierr));
 53: }

 55: static PetscErrorCode ourtaoadmmregularizerconstraintjacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
 56: {
 57:   PetscObjectUseFortranCallback(tao, _cb.regjacobian, (Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), (&tao, &x, &J, &Jpre, _ctx, &ierr));
 58: }

 60: PETSC_EXTERN void taoadmmsetmisfitobjectiveandgradientroutine_(Tao *tao, void (*func)(Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 61: {
 62:   CHKFORTRANNULLFUNCTION(func);
 63:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitobjgrad, (PetscVoidFn *)func, ctx);
 64:   if (!*ierr) *ierr = TaoADMMSetMisfitObjectiveAndGradientRoutine(*tao, ourtaoadmmmisfitobjgradroutine, ctx);
 65: }

 67: PETSC_EXTERN void taoadmmsetmisfithessianroutine_(Tao *tao, Mat *H, Mat *Hpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 68: {
 69:   CHKFORTRANNULLFUNCTION(func);
 70:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfithess, (PetscVoidFn *)func, ctx);
 71:   if (!*ierr) *ierr = TaoADMMSetMisfitHessianRoutine(*tao, *H, *Hpre, ourtaoadmmmisfithessroutine, ctx);
 72: }

 74: PETSC_EXTERN void taoadmmsetmisfitconstraintjacobian_(Tao *tao, Mat *J, Mat *Jpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 75: {
 76:   CHKFORTRANNULLFUNCTION(func);
 77:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitjacobian, (PetscVoidFn *)func, ctx);
 78:   if (!*ierr) *ierr = TaoADMMSetMisfitConstraintJacobian(*tao, *J, *Jpre, ourtaoadmmmisfitconstraintjacobian, ctx);
 79: }

 81: PETSC_EXTERN void taoadmmsetregularizerobjectiveandgradientroutine_(Tao *tao, void (*func)(Tao *, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 82: {
 83:   CHKFORTRANNULLFUNCTION(func);
 84:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.regobjgrad, (PetscVoidFn *)func, ctx);
 85:   if (!*ierr) *ierr = TaoADMMSetRegularizerObjectiveAndGradientRoutine(*tao, ourtaoadmmregularizerobjgradroutine, ctx);
 86: }

 88: PETSC_EXTERN void taoadmmsetregularizerhessianroutine_(Tao *tao, Mat *H, Mat *Hpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 89: {
 90:   CHKFORTRANNULLFUNCTION(func);
 91:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.reghess, (PetscVoidFn *)func, ctx);
 92:   if (!*ierr) *ierr = TaoADMMSetRegularizerHessianRoutine(*tao, *H, *Hpre, ourtaoadmmregularizerhessroutine, ctx);
 93: }

 95: PETSC_EXTERN void taoadmmsetregularizerconstraintjacobian_(Tao *tao, Mat *J, Mat *Jpre, void (*func)(Tao *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
 96: {
 97:   CHKFORTRANNULLFUNCTION(func);
 98:   *ierr = PetscObjectSetFortranCallback((PetscObject)*tao, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.misfitjacobian, (PetscVoidFn *)func, ctx);
 99:   if (!*ierr) *ierr = TaoADMMSetRegularizerConstraintJacobian(*tao, *J, *Jpre, ourtaoadmmregularizerconstraintjacobian, ctx);
100: }