1 #include <petsc/private/fortranimpl.h> 2 #include <petsc/private/taoimpl.h> 3 4 5 #if defined(PETSC_HAVE_FORTRAN_CAPS) 6 #define taosetobjectiveroutine_ TAOSETOBJECTIVEROUTINE 7 #define taosetgradientroutine_ TAOSETGRADIENTROUTINE 8 #define taosetobjectiveandgradientroutine_ TAOSETOBJECTIVEANDGRADIENTROUTINE 9 #define taosethessianroutine_ TAOSETHESSIANROUTINE 10 #define taosetseparableobjectiveroutine_ TAOSETSEPARABLEOBJECTIVEROUTINE 11 #define taosetjacobianroutine_ TAOSETJACOBIANROUTINE 12 #define taosetjacobianstateroutine_ TAOSETJACOBIANSTATEROUTINE 13 #define taosetjacobiandesignroutine_ TAOSETJACOBIANDESIGNROUTINE 14 #define taosetjacobianinequalityroutine_ TAOSETJACOBIANINEQUALITYROUTINE 15 #define taosetjacobianequalityroutine_ TAOSETJACOBIANEQUALITYROUTINE 16 #define taosetinequalityconstraintsroutine_ TAOSETINEQUALITYCONSTRAINTSROUTINE 17 #define taosetequalityconstraintsroutine_ TAOSETEQUALITYCONSTRAINTSROUTINE 18 #define taosetvariableboundsroutine_ TAOSETVARIABLEBOUNDSROUTINE 19 #define taosetconstraintsroutine_ TAOSETCONSTRAINTSROUTINE 20 #define taosetmonitor_ TAOSETMONITOR 21 #define taosettype_ TAOSETTYPE 22 #define taoview_ TAOVIEW 23 #define taogetconvergencehistory_ TAOGETCONVERGENCEHISTORY 24 #define taosetconvergencetest_ TAOSETCONVERGENCETEST 25 #define taogetoptionsprefix_ TAOGETOPTIONSPREFIX 26 #define taosetoptionsprefix_ TAOSETOPTIONSPREFIX 27 #define taoappendoptionsprefix_ TAOAPPENDOPTIONSPREFIX 28 #define taogettype_ TAOGETTYPE 29 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 30 31 #define taosetobjectiveroutine_ taosetobjectiveroutine 32 #define taosetgradientroutine_ taosetgradientroutine 33 #define taosetobjectiveandgradientroutine_ taosetobjectiveandgradientroutine 34 #define taosethessianroutine_ taosethessianroutine 35 #define taosetseparableobjectiveroutine_ taosetseparableobjectiveroutine 36 #define taosetjacobianroutine_ taosetjacobianroutine 37 #define taosetjacobianstateroutine_ taosetjacobianstateroutine 38 #define taosetjacobiandesignroutine_ taosetjacobiandesignroutine 39 #define taosetjacobianinequalityroutine_ taosetjacobianinequalityroutine 40 #define taosetjacobianequalityroutine_ taosetjacobianequalityroutine 41 #define taosetinequalityconstraintsroutine_ taosetinequalityconstraintsroutine 42 #define taosetequalityconstraintsroutine_ taosetequalityconstraintsroutine 43 #define taosetvariableboundsroutine_ taosetvariableboundsroutine 44 #define taosetconstraintsroutine_ taosetconstraintsroutine 45 #define taosetmonitor_ taosetmonitor 46 #define taosettype_ taosettype 47 #define taoview_ taoview 48 #define taogetconvergencehistory_ taogetconvergencehistory 49 #define taosetconvergencetest_ taosetconvergencetest 50 #define taogetoptionsprefix_ taogetoptionsprefix 51 #define taosetoptionsprefix_ taosetoptionsprefix 52 #define taoappendoptionsprefix_ taoappendoptionsprefix 53 #define taogettype_ taogettype 54 #endif 55 56 static int OBJ=0; /* objective routine index */ 57 static int GRAD=1; /* gradient routine index */ 58 static int OBJGRAD=2; /* objective and gradient routine */ 59 static int HESS=3; /* hessian routine index */ 60 static int SEPOBJ=4; /* separable objective routine index */ 61 static int JAC=5; /* jacobian routine index */ 62 static int JACSTATE=6; /* jacobian state routine index */ 63 static int JACDESIGN=7; /* jacobian design routine index */ 64 static int BOUNDS=8; 65 static int MON=9; /* monitor routine index */ 66 static int MONCTX=10; /* monitor routine index */ 67 static int MONDESTROY=11; /* monitor destroy index */ 68 static int CONVTEST=12; /* */ 69 static int CONSTRAINTS=13; 70 static int JACINEQ=14; 71 static int JACEQ=15; 72 static int CONINEQ=16; 73 static int CONEQ=17; 74 static int NFUNCS=18; 75 76 static PetscErrorCode ourtaoobjectiveroutine(Tao tao, Vec x, PetscReal *f, void *ctx) 77 { 78 PetscErrorCode ierr = 0; 79 (*(void (PETSC_STDCALL *)(Tao*,Vec*,PetscReal*,void*,PetscErrorCode*)) 80 (((PetscObject)tao)->fortran_func_pointers[OBJ]))(&tao,&x,f,ctx,&ierr); 81 CHKERRQ(ierr); 82 return 0; 83 } 84 85 static PetscErrorCode ourtaogradientroutine(Tao tao, Vec x, Vec g, void *ctx) 86 { 87 PetscErrorCode ierr = 0; 88 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 89 (((PetscObject)tao)->fortran_func_pointers[GRAD]))(&tao,&x,&g,ctx,&ierr); 90 CHKERRQ(ierr); 91 return 0; 92 93 } 94 95 static PetscErrorCode ourtaoobjectiveandgradientroutine(Tao tao, Vec x, PetscReal *f, Vec g, void* ctx) 96 { 97 PetscErrorCode ierr = 0; 98 (*(void (PETSC_STDCALL *)(Tao*,Vec*,PetscReal*,Vec*,void*,PetscErrorCode*)) 99 (((PetscObject)tao)->fortran_func_pointers[OBJGRAD]))(&tao,&x,f,&g,ctx,&ierr); 100 CHKERRQ(ierr); 101 return 0; 102 } 103 104 static PetscErrorCode ourtaohessianroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx) 105 { 106 PetscErrorCode ierr = 0; 107 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,Mat*,void*,PetscErrorCode*)) 108 (((PetscObject)tao)->fortran_func_pointers[HESS]))(&tao,&x,&H,&Hpre,ctx,&ierr); CHKERRQ(ierr); 109 return 0; 110 } 111 112 static PetscErrorCode ourtaojacobianroutine(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx) 113 { 114 PetscErrorCode ierr = 0; 115 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,Mat*,void*,PetscErrorCode*)) 116 (((PetscObject)tao)->fortran_func_pointers[JAC]))(&tao,&x,&H,&Hpre,ctx,&ierr); CHKERRQ(ierr); 117 return 0; 118 } 119 120 static PetscErrorCode ourtaojacobianstateroutine(Tao tao, Vec x, Mat H, Mat Hpre, Mat Hinv, void *ctx) 121 { 122 PetscErrorCode ierr = 0; 123 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,Mat*,Mat*,void*,PetscErrorCode*)) 124 (((PetscObject)tao)->fortran_func_pointers[JACSTATE]))(&tao,&x,&H,&Hpre,&Hinv,ctx,&ierr); CHKERRQ(ierr); 125 return 0; 126 } 127 128 static PetscErrorCode ourtaojacobiandesignroutine(Tao tao, Vec x, Mat H, void *ctx) 129 { 130 PetscErrorCode ierr = 0; 131 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,void*,PetscErrorCode*)) 132 (((PetscObject)tao)->fortran_func_pointers[JACDESIGN]))(&tao,&x,&H,ctx,&ierr); CHKERRQ(ierr); 133 return 0; 134 } 135 136 static PetscErrorCode ourtaoboundsroutine(Tao tao, Vec xl, Vec xu, void *ctx) 137 { 138 PetscErrorCode ierr = 0; 139 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 140 (((PetscObject)tao)->fortran_func_pointers[BOUNDS]))(&tao,&xl,&xu,ctx,&ierr); CHKERRQ(ierr); 141 return 0; 142 } 143 static PetscErrorCode ourtaoseparableobjectiveroutine(Tao tao, Vec x, Vec f, void *ctx) 144 { 145 PetscErrorCode ierr = 0; 146 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 147 (((PetscObject)tao)->fortran_func_pointers[SEPOBJ]))(&tao,&x,&f,ctx,&ierr); 148 return 0; 149 } 150 151 static PetscErrorCode ourtaomonitor(Tao tao, void *ctx) 152 { 153 PetscErrorCode ierr = 0; 154 (*(void (PETSC_STDCALL *)(Tao *, void*, PetscErrorCode*)) 155 (((PetscObject)tao)->fortran_func_pointers[MON]))(&tao,ctx,&ierr); 156 CHKERRQ(ierr); 157 return 0; 158 } 159 160 static PetscErrorCode ourtaomondestroy(void **ctx) 161 { 162 PetscErrorCode ierr = 0; 163 Tao tao = *(Tao*)ctx; 164 (*(void (PETSC_STDCALL *)(void*,PetscErrorCode*))(((PetscObject)tao)->fortran_func_pointers[MONDESTROY])) 165 ((void*)(PETSC_UINTPTR_T)((PetscObject)tao)->fortran_func_pointers[MONCTX],&ierr); 166 CHKERRQ(ierr); 167 return 0; 168 } 169 static PetscErrorCode ourtaoconvergencetest(Tao tao, void *ctx) 170 { 171 PetscErrorCode ierr = 0; 172 (*(void (PETSC_STDCALL *)(Tao *, void*, PetscErrorCode*)) 173 (((PetscObject)tao)->fortran_func_pointers[CONVTEST]))(&tao,ctx,&ierr); 174 CHKERRQ(ierr); 175 return 0; 176 } 177 178 179 static PetscErrorCode ourtaoconstraintsroutine(Tao tao, Vec x, Vec c, void *ctx) 180 { 181 PetscErrorCode ierr = 0; 182 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 183 (((PetscObject)tao)->fortran_func_pointers[CONSTRAINTS]))(&tao,&x,&c,ctx,&ierr); 184 CHKERRQ(ierr); 185 return 0; 186 187 } 188 189 static PetscErrorCode ourtaojacobianinequalityroutine(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx) 190 { 191 PetscErrorCode ierr = 0; 192 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,Mat*,void*,PetscErrorCode*)) 193 (((PetscObject)tao)->fortran_func_pointers[JACINEQ]))(&tao,&x,&J,&Jpre,ctx,&ierr); CHKERRQ(ierr); 194 return 0; 195 } 196 197 static PetscErrorCode ourtaojacobianequalityroutine(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx) 198 { 199 PetscErrorCode ierr = 0; 200 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Mat*,Mat*,void*,PetscErrorCode*)) 201 (((PetscObject)tao)->fortran_func_pointers[JACEQ]))(&tao,&x,&J,&Jpre,ctx,&ierr); CHKERRQ(ierr); 202 return 0; 203 } 204 205 static PetscErrorCode ourtaoinequalityconstraintsroutine(Tao tao, Vec x, Vec c, void *ctx) 206 { 207 PetscErrorCode ierr = 0; 208 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 209 (((PetscObject)tao)->fortran_func_pointers[CONINEQ]))(&tao,&x,&c,ctx,&ierr); 210 CHKERRQ(ierr); 211 return 0; 212 213 } 214 215 static PetscErrorCode ourtaoequalityconstraintsroutine(Tao tao, Vec x, Vec c, void *ctx) 216 { 217 PetscErrorCode ierr = 0; 218 (*(void (PETSC_STDCALL *)(Tao*,Vec*,Vec*,void*,PetscErrorCode*)) 219 (((PetscObject)tao)->fortran_func_pointers[CONEQ]))(&tao,&x,&c,ctx,&ierr); 220 CHKERRQ(ierr); 221 return 0; 222 223 } 224 225 226 EXTERN_C_BEGIN 227 228 229 void PETSC_STDCALL taosetobjectiveroutine_(Tao *tao, void (PETSC_STDCALL *func)(Tao*, Vec *, PetscReal *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 230 { 231 CHKFORTRANNULLOBJECT(ctx); 232 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 233 if (!func) { 234 *ierr = TaoSetObjectiveRoutine(*tao,0,ctx); 235 } else { 236 ((PetscObject)*tao)->fortran_func_pointers[OBJ] = (PetscVoidFunction)func; 237 *ierr = TaoSetObjectiveRoutine(*tao, ourtaoobjectiveroutine,ctx); 238 } 239 } 240 241 void PETSC_STDCALL taosetgradientroutine_(Tao *tao, void (PETSC_STDCALL *func)(Tao*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 242 { 243 CHKFORTRANNULLOBJECT(ctx); 244 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 245 if (!func) { 246 *ierr = TaoSetGradientRoutine(*tao,0,ctx); 247 } else { 248 ((PetscObject)*tao)->fortran_func_pointers[GRAD] = (PetscVoidFunction)func; 249 *ierr = TaoSetGradientRoutine(*tao, ourtaogradientroutine,ctx); 250 } 251 } 252 253 void PETSC_STDCALL taosetobjectiveandgradientroutine_(Tao *tao, void (PETSC_STDCALL *func)(Tao*, Vec *, PetscReal *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 254 { 255 CHKFORTRANNULLOBJECT(ctx); 256 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 257 if (!func) { 258 *ierr = TaoSetObjectiveAndGradientRoutine(*tao,0,ctx); 259 } else { 260 ((PetscObject)*tao)->fortran_func_pointers[OBJGRAD] = (PetscVoidFunction)func; 261 *ierr = TaoSetObjectiveAndGradientRoutine(*tao, ourtaoobjectiveandgradientroutine,ctx); 262 } 263 } 264 265 266 267 268 void PETSC_STDCALL taosetseparableobjectiveroutine_(Tao *tao, Vec *F, void (PETSC_STDCALL *func)(Tao*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 269 { 270 CHKFORTRANNULLOBJECT(ctx); 271 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 272 if (!func) { 273 *ierr = TaoSetSeparableObjectiveRoutine(*tao,*F,0,ctx); 274 } else { 275 ((PetscObject)*tao)->fortran_func_pointers[SEPOBJ] = (PetscVoidFunction)func; 276 *ierr = TaoSetSeparableObjectiveRoutine(*tao,*F, ourtaoseparableobjectiveroutine,ctx); 277 } 278 } 279 280 281 282 void PETSC_STDCALL taosetjacobianroutine_(Tao *tao, Mat *J, Mat *Jp, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 283 { 284 CHKFORTRANNULLOBJECT(ctx); 285 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 286 if (!func) { 287 *ierr = TaoSetJacobianRoutine(*tao,*J,*Jp,0,ctx); 288 } else { 289 ((PetscObject)*tao)->fortran_func_pointers[JAC] = (PetscVoidFunction)func; 290 *ierr = TaoSetJacobianRoutine(*tao,*J, *Jp, ourtaojacobianroutine,ctx); 291 } 292 } 293 294 void PETSC_STDCALL taosetjacobianstateroutine_(Tao *tao, Mat *J, Mat *Jp, Mat*Jinv, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, Mat *, Mat*, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 295 { 296 CHKFORTRANNULLOBJECT(ctx); 297 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 298 if (!func) { 299 *ierr = TaoSetJacobianStateRoutine(*tao,*J,*Jp,*Jinv,0,ctx); 300 } else { 301 ((PetscObject)*tao)->fortran_func_pointers[JACSTATE] = (PetscVoidFunction)func; 302 *ierr = TaoSetJacobianStateRoutine(*tao,*J, *Jp, *Jinv, ourtaojacobianstateroutine,ctx); 303 } 304 } 305 306 void PETSC_STDCALL taosetjacobiandesignroutine_(Tao *tao, Mat *J, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 307 { 308 CHKFORTRANNULLOBJECT(ctx); 309 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 310 if (!func) { 311 *ierr = TaoSetJacobianDesignRoutine(*tao,*J,0,ctx); 312 } else { 313 ((PetscObject)*tao)->fortran_func_pointers[JACDESIGN] = (PetscVoidFunction)func; 314 *ierr = TaoSetJacobianDesignRoutine(*tao,*J, ourtaojacobiandesignroutine,ctx); 315 } 316 } 317 318 319 void PETSC_STDCALL taosethessianroutine_(Tao *tao, Mat *J, Mat *Jp, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, Mat *,void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 320 { 321 CHKFORTRANNULLOBJECT(ctx); 322 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 323 if (!func) { 324 *ierr = TaoSetHessianRoutine(*tao,*J,*Jp,0,ctx); 325 } else { 326 ((PetscObject)*tao)->fortran_func_pointers[HESS] = (PetscVoidFunction)func; 327 *ierr = TaoSetHessianRoutine(*tao,*J, *Jp, ourtaohessianroutine,ctx); 328 } 329 } 330 331 void PETSC_STDCALL taosetvariableboundsroutine_(Tao *tao, void (PETSC_STDCALL *func)(Tao*,Vec*,Vec*,void*,PetscErrorCode*),void *ctx, PetscErrorCode *ierr) 332 { 333 CHKFORTRANNULLOBJECT(ctx); 334 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 335 if (func) { 336 ((PetscObject)*tao)->fortran_func_pointers[BOUNDS] = (PetscVoidFunction)func; 337 *ierr = TaoSetVariableBoundsRoutine(*tao,ourtaoboundsroutine,ctx); 338 } else { 339 *ierr = TaoSetVariableBoundsRoutine(*tao,0,ctx); 340 } 341 342 } 343 void PETSC_STDCALL taosetmonitor_(Tao *tao, void (PETSC_STDCALL *func)(Tao*,void*,PetscErrorCode*),void *ctx, void (PETSC_STDCALL *mondestroy)(void*,PetscErrorCode*),PetscErrorCode *ierr) 344 { 345 CHKFORTRANNULLOBJECT(ctx); 346 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 347 if (func) { 348 ((PetscObject)*tao)->fortran_func_pointers[MON] = (PetscVoidFunction)func; 349 if (FORTRANNULLFUNCTION(mondestroy)){ 350 *ierr = TaoSetMonitor(*tao,ourtaomonitor,*tao,NULL); 351 } else { 352 *ierr = TaoSetMonitor(*tao,ourtaomonitor,*tao,ourtaomondestroy); 353 } 354 } 355 } 356 357 void PETSC_STDCALL taosetconvergencetest_(Tao *tao, void (PETSC_STDCALL *func)(Tao*,void*,PetscErrorCode*),void *ctx, PetscErrorCode *ierr) 358 { 359 CHKFORTRANNULLOBJECT(ctx); 360 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 361 if (!func) { 362 *ierr = TaoSetConvergenceTest(*tao,0,ctx); 363 } else { 364 ((PetscObject)*tao)->fortran_func_pointers[CONVTEST] = (PetscVoidFunction)func; 365 *ierr = TaoSetConvergenceTest(*tao,ourtaoconvergencetest,ctx); 366 } 367 } 368 369 370 void PETSC_STDCALL taosetconstraintsroutine_(Tao *tao, Vec *C, void (PETSC_STDCALL *func)(Tao*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 371 { 372 CHKFORTRANNULLOBJECT(ctx); 373 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 374 if (!func) { 375 *ierr = TaoSetConstraintsRoutine(*tao,*C,0,ctx); 376 } else { 377 ((PetscObject)*tao)->fortran_func_pointers[CONSTRAINTS] = (PetscVoidFunction)func; 378 *ierr = TaoSetConstraintsRoutine(*tao, *C, ourtaoconstraintsroutine,ctx); 379 } 380 } 381 382 383 void PETSC_STDCALL taosettype_(Tao *tao, CHAR type_name PETSC_MIXED_LEN(len), PetscErrorCode *ierr PETSC_END_LEN(len)) 384 385 { 386 char *t; 387 388 FIXCHAR(type_name,len,t); 389 *ierr = TaoSetType(*tao,t); 390 FREECHAR(type_name,t); 391 392 } 393 394 void PETSC_STDCALL taoview_(Tao *tao, PetscViewer *viewer, PetscErrorCode *ierr) 395 { 396 PetscViewer v; 397 PetscPatchDefaultViewers_Fortran(viewer,v); 398 *ierr = TaoView(*tao,v); 399 } 400 401 void PETSC_STDCALL taogetconvergencehistory_(Tao *tao, PetscInt *nhist, PetscErrorCode *ierr) 402 { 403 *ierr = TaoGetConvergenceHistory(*tao,NULL,NULL,NULL,NULL,nhist); 404 } 405 406 void PETSC_STDCALL taogetoptionsprefix_(Tao *tao, CHAR prefix PETSC_MIXED_LEN(len), PetscErrorCode *ierr PETSC_END_LEN(len)) 407 { 408 const char *name; 409 *ierr = TaoGetOptionsPrefix(*tao,&name); 410 *ierr = PetscStrncpy(prefix,name,len); if (*ierr) return; 411 FIXRETURNCHAR(PETSC_TRUE,prefix,len); 412 413 } 414 415 void PETSC_STDCALL taoappendoptionsprefix_(Tao *tao, CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 416 { 417 char *name; 418 FIXCHAR(prefix,len,name); 419 *ierr = TaoAppendOptionsPrefix(*tao,name); 420 FREECHAR(prefix,name); 421 } 422 423 void PETSC_STDCALL taosetoptionsprefix_(Tao *tao, CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 424 { 425 char *t; 426 FIXCHAR(prefix,len,t); 427 *ierr = TaoSetOptionsPrefix(*tao,t); 428 FREECHAR(prefix,t); 429 } 430 431 void PETSC_STDCALL taogettype_(Tao *tao, CHAR name PETSC_MIXED_LEN(len), PetscErrorCode *ierr PETSC_END_LEN(len)) 432 { 433 const char *tname; 434 *ierr = TaoGetType(*tao,&tname); 435 *ierr = PetscStrncpy(name,tname,len); if (*ierr) return; 436 FIXRETURNCHAR(PETSC_TRUE,name,len); 437 438 } 439 440 441 void PETSC_STDCALL taosetjacobianinequalityroutine_(Tao *tao, Mat *J, Mat *Jp, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, Mat *,void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 442 { 443 CHKFORTRANNULLOBJECT(ctx); 444 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 445 if (!func) { 446 *ierr = TaoSetJacobianInequalityRoutine(*tao,*J,*Jp,0,ctx); 447 } else { 448 ((PetscObject)*tao)->fortran_func_pointers[JACINEQ] = (PetscVoidFunction)func; 449 *ierr = TaoSetJacobianInequalityRoutine(*tao,*J, *Jp, ourtaojacobianinequalityroutine,ctx); 450 } 451 } 452 453 void PETSC_STDCALL taosetjacobianequalityroutine_(Tao *tao, Mat *J, Mat *Jp, void (PETSC_STDCALL *func)(Tao*, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 454 { 455 CHKFORTRANNULLOBJECT(ctx); 456 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 457 if (!func) { 458 *ierr = TaoSetJacobianEqualityRoutine(*tao,*J,*Jp,0,ctx); 459 } else { 460 ((PetscObject)*tao)->fortran_func_pointers[JACEQ] = (PetscVoidFunction)func; 461 *ierr = TaoSetJacobianEqualityRoutine(*tao,*J, *Jp, ourtaojacobianequalityroutine,ctx); 462 } 463 } 464 465 466 void PETSC_STDCALL taosetinequalityconstraintsroutine_(Tao *tao, Vec *C, void (PETSC_STDCALL *func)(Tao*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 467 { 468 CHKFORTRANNULLOBJECT(ctx); 469 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 470 if (!func) { 471 *ierr = TaoSetInequalityConstraintsRoutine(*tao,*C,0,ctx); 472 } else { 473 ((PetscObject)*tao)->fortran_func_pointers[CONINEQ] = (PetscVoidFunction)func; 474 *ierr = TaoSetInequalityConstraintsRoutine(*tao, *C, ourtaoinequalityconstraintsroutine,ctx); 475 } 476 } 477 478 void PETSC_STDCALL taosetequalityconstraintsroutine_(Tao *tao, Vec *C, void (PETSC_STDCALL *func)(Tao*, Vec *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 479 { 480 CHKFORTRANNULLOBJECT(ctx); 481 PetscObjectAllocateFortranPointers(*tao,NFUNCS); 482 if (!func) { 483 *ierr = TaoSetEqualityConstraintsRoutine(*tao,*C,0,ctx); 484 } else { 485 ((PetscObject)*tao)->fortran_func_pointers[CONEQ] = (PetscVoidFunction)func; 486 *ierr = TaoSetEqualityConstraintsRoutine(*tao, *C, ourtaoequalityconstraintsroutine,ctx); 487 } 488 } 489 490 491 EXTERN_C_END 492 493 494