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 taolinesearchsettype_ TAOLINESEARCHSETTYPE 11 #define taolinesearchviewfromoptions_ TAOLINESEARCHVIEWFROMOPTIONS 12 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 13 14 #define taolinesearchsetobjectiveroutine_ taolinesearchsetobjectiveroutine 15 #define taolinesearchsetgradientroutine_ taolinesearchsetgradientroutine 16 #define taolinesearchsetobjectiveandgradientroutine_ taolinesearchsetobjectiveandgradientroutine 17 #define taolinesearchsetobjectiveandgtsroutine_ taolinesearchsetobjectiveandgtsroutine 18 #define taolinesearchview_ taolinesearchview 19 #define taolinesearchsettype_ taolinesearchsettype 20 #define taolinesearchviewfromoptions_ taolinesearchviewfromoptions 21 #endif 22 23 static int OBJ=0; 24 static int GRAD=1; 25 static int OBJGRAD=2; 26 static int OBJGTS=3; 27 static int NFUNCS=4; 28 29 static PetscErrorCode ourtaolinesearchobjectiveroutine(TaoLineSearch ls, Vec x, PetscReal *f, void *ctx) 30 { 31 PetscErrorCode ierr = 0; 32 (*(void (*)(TaoLineSearch*,Vec*,PetscReal*,void*,PetscErrorCode*)) (((PetscObject)ls)->fortran_func_pointers[OBJ]))(&ls,&x,f,ctx,&ierr);CHKERRQ(ierr); 33 return 0; 34 } 35 36 static PetscErrorCode ourtaolinesearchgradientroutine(TaoLineSearch ls, Vec x, Vec g, void *ctx) 37 { 38 PetscErrorCode ierr = 0; 39 (*(void (*)(TaoLineSearch*,Vec*,Vec*,void*,PetscErrorCode*)) (((PetscObject)ls)->fortran_func_pointers[GRAD]))(&ls,&x,&g,ctx,&ierr);CHKERRQ(ierr); 40 return 0; 41 42 } 43 44 static PetscErrorCode ourtaolinesearchobjectiveandgradientroutine(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, void* ctx) 45 { 46 PetscErrorCode ierr = 0; 47 (*(void (*)(TaoLineSearch*,Vec*,PetscReal*,Vec*,void*,PetscErrorCode*))(((PetscObject)ls)->fortran_func_pointers[OBJGRAD]))(&ls,&x,f,&g,ctx,&ierr);CHKERRQ(ierr); 48 return 0; 49 } 50 51 static PetscErrorCode ourtaolinesearchobjectiveandgtsroutine(TaoLineSearch ls, Vec x, Vec s, PetscReal *f, PetscReal *gts, void* ctx) 52 { 53 PetscErrorCode ierr = 0; 54 (*(void (*)(TaoLineSearch*,Vec*,Vec*,PetscReal*,PetscReal*,void*,PetscErrorCode*)) (((PetscObject)ls)->fortran_func_pointers[OBJGTS]))(&ls,&x,&s,f,gts,ctx,&ierr);CHKERRQ(ierr); 55 return 0; 56 } 57 58 PETSC_EXTERN void taolinesearchsetobjectiveroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch*, Vec *, PetscReal *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 59 { 60 PetscObjectAllocateFortranPointers(*ls,NFUNCS); 61 if (!func) { 62 *ierr = TaoLineSearchSetObjectiveRoutine(*ls,0,ctx); 63 } else { 64 ((PetscObject)*ls)->fortran_func_pointers[OBJ] = (PetscVoidFunction)func; 65 *ierr = TaoLineSearchSetObjectiveRoutine(*ls, ourtaolinesearchobjectiveroutine,ctx); 66 } 67 } 68 69 PETSC_EXTERN void taolinesearchsetgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 70 { 71 PetscObjectAllocateFortranPointers(*ls,NFUNCS); 72 if (!func) { 73 *ierr = TaoLineSearchSetGradientRoutine(*ls,0,ctx); 74 } else { 75 ((PetscObject)*ls)->fortran_func_pointers[GRAD] = (PetscVoidFunction)func; 76 *ierr = TaoLineSearchSetGradientRoutine(*ls, ourtaolinesearchgradientroutine,ctx); 77 } 78 } 79 80 PETSC_EXTERN void taolinesearchsetobjectiveandgradientroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch*, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 81 { 82 PetscObjectAllocateFortranPointers(*ls,NFUNCS); 83 if (!func) { 84 *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls,0,ctx); 85 } else { 86 ((PetscObject)*ls)->fortran_func_pointers[OBJGRAD] = (PetscVoidFunction)func; 87 *ierr = TaoLineSearchSetObjectiveAndGradientRoutine(*ls, ourtaolinesearchobjectiveandgradientroutine,ctx); 88 } 89 } 90 91 PETSC_EXTERN void taolinesearchsetobjectiveandgtsroutine_(TaoLineSearch *ls, void (*func)(TaoLineSearch*, Vec *, Vec *, PetscReal*, PetscReal*,void*, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 92 { 93 PetscObjectAllocateFortranPointers(*ls,NFUNCS); 94 if (!func) { 95 *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls,0,ctx); 96 } else { 97 ((PetscObject)*ls)->fortran_func_pointers[OBJGTS] = (PetscVoidFunction)func; 98 *ierr = TaoLineSearchSetObjectiveAndGTSRoutine(*ls, ourtaolinesearchobjectiveandgtsroutine,ctx); 99 } 100 } 101 102 PETSC_EXTERN void taolinesearchsettype_(TaoLineSearch *ls, char* type_name, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 103 104 { 105 char *t; 106 107 FIXCHAR(type_name,len,t); 108 *ierr = TaoLineSearchSetType(*ls,t);if (*ierr) return; 109 FREECHAR(type_name,t); 110 111 } 112 113 PETSC_EXTERN void taolinesearchview_(TaoLineSearch *ls, PetscViewer *viewer, PetscErrorCode *ierr) 114 { 115 PetscViewer v; 116 PetscPatchDefaultViewers_Fortran(viewer,v); 117 *ierr = TaoLineSearchView(*ls,v); 118 } 119 120 PETSC_EXTERN void taolinesearchgetoptionsprefix_(TaoLineSearch *ls, char* prefix, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 121 { 122 const char *name; 123 *ierr = TaoLineSearchGetOptionsPrefix(*ls,&name); 124 *ierr = PetscStrncpy(prefix,name,len); if (*ierr) return; 125 FIXRETURNCHAR(PETSC_TRUE,prefix,len); 126 127 } 128 129 PETSC_EXTERN void taolinesearchappendoptionsprefix_(TaoLineSearch *ls, char* prefix,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 130 { 131 char *name; 132 FIXCHAR(prefix,len,name); 133 *ierr = TaoLineSearchAppendOptionsPrefix(*ls,name);if (*ierr) return; 134 FREECHAR(prefix,name); 135 } 136 137 PETSC_EXTERN void taolinesearchsetoptionsprefix_(TaoLineSearch *ls, char* prefix,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 138 { 139 char *t; 140 FIXCHAR(prefix,len,t); 141 *ierr = TaoLineSearchSetOptionsPrefix(*ls,t);if (*ierr) return; 142 FREECHAR(prefix,t); 143 } 144 145 PETSC_EXTERN void taolinesearchgettype_(TaoLineSearch *ls, char* name, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 146 { 147 const char *tname; 148 *ierr = TaoLineSearchGetType(*ls,&tname); 149 *ierr = PetscStrncpy(name,tname,len); if (*ierr) return; 150 FIXRETURNCHAR(PETSC_TRUE,name,len); 151 152 } 153 PETSC_EXTERN void taolinesearchviewfromoptions_(TaoLineSearch *ao,PetscObject obj,char* type,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 154 { 155 char *t; 156 157 FIXCHAR(type,len,t); 158 CHKFORTRANNULLOBJECT(obj); 159 *ierr = TaoLineSearchViewFromOptions(*ao,obj,t);if (*ierr) return; 160 FREECHAR(type,t); 161 } 162