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