1cdef extern from * nogil: 2 3 ctypedef const char* PetscTAOType "TaoType" 4 PetscTAOType TAOLMVM 5 PetscTAOType TAONLS 6 PetscTAOType TAONTR 7 PetscTAOType TAONTL 8 PetscTAOType TAOCG 9 PetscTAOType TAOTRON 10 PetscTAOType TAOOWLQN 11 PetscTAOType TAOBMRM 12 PetscTAOType TAOBLMVM 13 PetscTAOType TAOBQNLS 14 PetscTAOType TAOBNCG 15 PetscTAOType TAOBNLS 16 PetscTAOType TAOBNTR 17 PetscTAOType TAOBNTL 18 PetscTAOType TAOBQNKLS 19 PetscTAOType TAOBQNKTR 20 PetscTAOType TAOBQNKTL 21 PetscTAOType TAOBQPIP 22 PetscTAOType TAOGPCG 23 PetscTAOType TAONM 24 PetscTAOType TAOPOUNDERS 25 PetscTAOType TAOBRGN 26 PetscTAOType TAOLCL 27 PetscTAOType TAOSSILS 28 PetscTAOType TAOSSFLS 29 PetscTAOType TAOASILS 30 PetscTAOType TAOASFLS 31 PetscTAOType TAOIPM 32 PetscTAOType TAOPDIPM 33 PetscTAOType TAOSHELL 34 PetscTAOType TAOADMM 35 PetscTAOType TAOALMM 36 PetscTAOType TAOPYTHON 37 38 ctypedef enum PetscTAOConvergedReason "TaoConvergedReason": 39 # iterating 40 TAO_CONTINUE_ITERATING 41 # converged 42 TAO_CONVERGED_GATOL 43 TAO_CONVERGED_GRTOL 44 TAO_CONVERGED_GTTOL 45 TAO_CONVERGED_STEPTOL 46 TAO_CONVERGED_MINF 47 TAO_CONVERGED_USER 48 # diverged 49 TAO_DIVERGED_MAXITS 50 TAO_DIVERGED_NAN 51 TAO_DIVERGED_MAXFCN 52 TAO_DIVERGED_LS_FAILURE 53 TAO_DIVERGED_TR_REDUCTION 54 TAO_DIVERGED_USER 55 56 ctypedef PetscErrorCode (*PetscTaoMonitorDestroy)(void**) 57 ctypedef PetscErrorCode PetscTaoConvergenceTest(PetscTAO, void*) except PETSC_ERR_PYTHON 58 ctypedef PetscErrorCode PetscTaoMonitor(PetscTAO, void*) except PETSC_ERR_PYTHON 59 ctypedef PetscErrorCode PetscTaoObjective(PetscTAO, PetscVec, PetscReal*, void*) except PETSC_ERR_PYTHON 60 ctypedef PetscErrorCode PetscTaoResidual(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 61 ctypedef PetscErrorCode PetscTaoGradient(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 62 ctypedef PetscErrorCode PetscTaoObjGrad(PetscTAO, PetscVec, PetscReal*, PetscVec, void*) except PETSC_ERR_PYTHON 63 ctypedef PetscErrorCode PetscTaoRegularizerObjGrad(PetscTAO, PetscVec, PetscReal*, PetscVec, void*) except PETSC_ERR_PYTHON 64 ctypedef PetscErrorCode PetscTaoVarBounds(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 65 ctypedef PetscErrorCode PetscTaoConstraints(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 66 ctypedef PetscErrorCode PetscTaoEqualityConstraints(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 67 ctypedef PetscErrorCode PetscTaoInequalityConstraints(PetscTAO, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 68 ctypedef PetscErrorCode PetscTaoHessian(PetscTAO, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 69 ctypedef PetscErrorCode PetscTaoRegularizerHessian(PetscTAO, PetscVec, PetscMat, void*) except PETSC_ERR_PYTHON 70 ctypedef PetscErrorCode PetscTaoJacobian(PetscTAO, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 71 ctypedef PetscErrorCode PetscTaoJacobianResidual(PetscTAO, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 72 ctypedef PetscErrorCode PetscTaoJacobianState(PetscTAO, PetscVec, PetscMat, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 73 ctypedef PetscErrorCode PetscTaoJacobianDesign(PetscTAO, PetscVec, PetscMat, void*) except PETSC_ERR_PYTHON 74 ctypedef PetscErrorCode PetscTaoJacobianEquality(PetscTAO, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 75 ctypedef PetscErrorCode PetscTaoJacobianInequality(PetscTAO, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 76 ctypedef PetscErrorCode PetscTaoUpdateFunction(PetscTAO, PetscInt, void*) except PETSC_ERR_PYTHON 77 78 ctypedef enum PetscTAOBNCGType "TaoBNCGType": 79 TAO_BNCG_GD 80 TAO_BNCG_PCGD 81 TAO_BNCG_HS 82 TAO_BNCG_FR 83 TAO_BNCG_PRP 84 TAO_BNCG_PRP_PLUS 85 TAO_BNCG_DY 86 TAO_BNCG_HZ 87 TAO_BNCG_DK 88 TAO_BNCG_KD 89 TAO_BNCG_SSML_BFGS 90 TAO_BNCG_SSML_DFP 91 TAO_BNCG_SSML_BRDN 92 93 ctypedef enum PetscTAOALMMType "TaoALMMType": 94 TAO_ALMM_CLASSIC 95 TAO_ALMM_PHR 96 97 PetscErrorCode TaoMonitor(PetscTAO, PetscInt, PetscReal, PetscReal, PetscReal, PetscReal) 98 PetscErrorCode TaoView(PetscTAO, PetscViewer) 99 PetscErrorCode TaoDestroy(PetscTAO*) 100 PetscErrorCode TaoCreate(MPI_Comm, PetscTAO*) 101 PetscErrorCode TaoSetOptionsPrefix(PetscTAO, char[]) 102 PetscErrorCode TaoAppendOptionsPrefix(PetscTAO, char[]) 103 PetscErrorCode TaoGetOptionsPrefix(PetscTAO, char*[]) 104 PetscErrorCode TaoSetFromOptions(PetscTAO) 105 PetscErrorCode TaoSetType(PetscTAO, PetscTAOType) 106 PetscErrorCode TaoGetType(PetscTAO, PetscTAOType*) 107 108 PetscErrorCode TaoSetUp(PetscTAO) 109 PetscErrorCode TaoSolve(PetscTAO) 110 111 PetscErrorCode TaoSetTolerances(PetscTAO, PetscReal, PetscReal, PetscReal) 112 PetscErrorCode TaoParametersInitialize(PetscTAO) 113 PetscErrorCode TaoGetTolerances(PetscTAO, PetscReal*, PetscReal*, PetscReal*) 114 PetscErrorCode TaoSetConstraintTolerances(PetscTAO, PetscReal, PetscReal) 115 PetscErrorCode TaoGetConstraintTolerances(PetscTAO, PetscReal*, PetscReal*) 116 117 PetscErrorCode TaoSetFunctionLowerBound(PetscTAO, PetscReal) 118 PetscErrorCode TaoSetMaximumIterations(PetscTAO, PetscInt) 119 PetscErrorCode TaoGetMaximumIterations(PetscTAO, PetscInt*) 120 PetscErrorCode TaoSetMaximumFunctionEvaluations(PetscTAO, PetscInt) 121 PetscErrorCode TaoGetMaximumFunctionEvaluations(PetscTAO, PetscInt*) 122 PetscErrorCode TaoSetIterationNumber(PetscTAO, PetscInt) 123 PetscErrorCode TaoGetIterationNumber(PetscTAO, PetscInt*) 124 125 PetscErrorCode TaoSetTrustRegionTolerance(PetscTAO, PetscReal) 126 PetscErrorCode TaoGetInitialTrustRegionRadius(PetscTAO, PetscReal*) 127 PetscErrorCode TaoGetTrustRegionRadius(PetscTAO, PetscReal*) 128 PetscErrorCode TaoSetTrustRegionRadius(PetscTAO, PetscReal) 129 130 PetscErrorCode TaoDefaultConvergenceTest(PetscTAO, void*) except PETSC_ERR_PYTHON 131 PetscErrorCode TaoSetConvergenceTest(PetscTAO, PetscTaoConvergenceTest*, void*) 132 PetscErrorCode TaoSetConvergedReason(PetscTAO, PetscTAOConvergedReason) 133 PetscErrorCode TaoGetConvergedReason(PetscTAO, PetscTAOConvergedReason*) 134 PetscErrorCode TaoLogConvergenceHistory(PetscTAO, PetscReal, PetscReal, PetscReal, PetscInt) 135 PetscErrorCode TaoGetSolutionStatus(PetscTAO, PetscInt*, 136 PetscReal*, PetscReal*, 137 PetscReal*, PetscReal*, 138 PetscTAOConvergedReason*) 139 140 PetscErrorCode TaoMonitorSet(PetscTAO, PetscTaoMonitor, void*, PetscTaoMonitorDestroy) 141 PetscErrorCode TaoMonitorCancel(PetscTAO) 142 143 PetscErrorCode TaoComputeObjective(PetscTAO, PetscVec, PetscReal*) 144 PetscErrorCode TaoComputeResidual(PetscTAO, PetscVec, PetscVec) 145 PetscErrorCode TaoComputeGradient(PetscTAO, PetscVec, PetscVec) 146 PetscErrorCode TaoComputeObjectiveAndGradient(PetscTAO, PetscVec, PetscReal*, PetscVec) 147 PetscErrorCode TaoComputeConstraints(PetscTAO, PetscVec, PetscVec) 148 PetscErrorCode TaoComputeDualVariables(PetscTAO, PetscVec, PetscVec) 149 PetscErrorCode TaoComputeVariableBounds(PetscTAO) 150 PetscErrorCode TaoComputeHessian(PetscTAO, PetscVec, PetscMat, PetscMat) 151 PetscErrorCode TaoComputeJacobian(PetscTAO, PetscVec, PetscMat, PetscMat) 152 153 PetscErrorCode TaoSetSolution(PetscTAO, PetscVec) 154 PetscErrorCode TaoSetConstraintsVec(PetscTAO, PetscVec) 155 PetscErrorCode TaoSetVariableBounds(PetscTAO, PetscVec, PetscVec) 156 157 PetscErrorCode TaoGetSolution(PetscTAO, PetscVec*) 158 PetscErrorCode TaoSetGradientNorm(PetscTAO, PetscMat) 159 PetscErrorCode TaoGetGradientNorm(PetscTAO, PetscMat*) 160 PetscErrorCode TaoLMVMSetH0(PetscTAO, PetscMat) 161 PetscErrorCode TaoLMVMGetH0(PetscTAO, PetscMat*) 162 PetscErrorCode TaoLMVMGetH0KSP(PetscTAO, PetscKSP*) 163 PetscErrorCode TaoBNCGGetType(PetscTAO, PetscTAOBNCGType*) 164 PetscErrorCode TaoBNCGSetType(PetscTAO, PetscTAOBNCGType) 165 PetscErrorCode TaoGetVariableBounds(PetscTAO, PetscVec*, PetscVec*) 166 167 PetscErrorCode TaoSetObjective(PetscTAO, PetscTaoObjective*, void*) 168 PetscErrorCode TaoSetGradient(PetscTAO, PetscVec, PetscTaoGradient*, void*) 169 PetscErrorCode TaoSetObjectiveAndGradient(PetscTAO, PetscVec, PetscTaoObjGrad*, void*) 170 PetscErrorCode TaoSetHessian(PetscTAO, PetscMat, PetscMat, PetscTaoHessian*, void*) 171 PetscErrorCode TaoGetObjective(PetscTAO, PetscTaoObjective**, void**) 172 PetscErrorCode TaoGetGradient(PetscTAO, PetscVec*, PetscTaoGradient**, void**) 173 PetscErrorCode TaoGetObjectiveAndGradient(PetscTAO, PetscVec*, PetscTaoObjGrad**, void**) 174 PetscErrorCode TaoGetHessian(PetscTAO, PetscMat*, PetscMat*, PetscTaoHessian**, void**) 175 PetscErrorCode TaoSetResidualRoutine(PetscTAO, PetscVec, PetscTaoResidual, void*) 176 PetscErrorCode TaoSetVariableBoundsRoutine(PetscTAO, PetscTaoVarBounds*, void*) 177 PetscErrorCode TaoSetConstraintsRoutine(PetscTAO, PetscVec, PetscTaoConstraints*, void*) 178 PetscErrorCode TaoSetJacobianRoutine(PetscTAO, PetscMat, PetscMat, PetscTaoJacobian*, void*) 179 PetscErrorCode TaoSetJacobianResidualRoutine(PetscTAO, PetscMat, PetscMat, PetscTaoJacobianResidual*, void*) 180 PetscErrorCode TaoSetStateDesignIS(PetscTAO, PetscIS, PetscIS) 181 PetscErrorCode TaoSetJacobianStateRoutine(PetscTAO, PetscMat, PetscMat, PetscMat, PetscTaoJacobianState*, void*) 182 PetscErrorCode TaoSetJacobianDesignRoutine(PetscTAO, PetscMat, PetscTaoJacobianDesign*, void*) 183 PetscErrorCode TaoGetLMVMMatrix(PetscTAO, PetscMat*) 184 PetscErrorCode TaoSetLMVMMatrix(PetscTAO, PetscMat) 185 186 PetscErrorCode TaoSetEqualityConstraintsRoutine(PetscTAO, PetscVec, PetscTaoEqualityConstraints*, void*) 187 PetscErrorCode TaoGetEqualityConstraintsRoutine(PetscTAO, PetscVec*, PetscTaoEqualityConstraints**, void**) 188 PetscErrorCode TaoSetJacobianEqualityRoutine(PetscTAO, PetscMat, PetscMat, PetscTaoJacobianEquality*, void*) 189 PetscErrorCode TaoGetJacobianEqualityRoutine(PetscTAO, PetscMat*, PetscMat*, PetscTaoJacobianEquality**, void**) 190 PetscErrorCode TaoSetInequalityConstraintsRoutine(PetscTAO, PetscVec, PetscTaoInequalityConstraints*, void*) 191 PetscErrorCode TaoGetInequalityConstraintsRoutine(PetscTAO, PetscVec*, PetscTaoInequalityConstraints**, void**) 192 PetscErrorCode TaoSetJacobianInequalityRoutine(PetscTAO, PetscMat, PetscMat, PetscTaoJacobianInequality*, void*) 193 PetscErrorCode TaoGetJacobianInequalityRoutine(PetscTAO, PetscMat*, PetscMat*, PetscTaoJacobianInequality**, void**) 194 PetscErrorCode TaoSetUpdate(PetscTAO, PetscTaoUpdateFunction*, void*) 195 196 PetscErrorCode TaoSetInitialTrustRegionRadius(PetscTAO, PetscReal) 197 198 PetscErrorCode TaoGetKSP(PetscTAO, PetscKSP*) 199 PetscErrorCode TaoGetLineSearch(PetscTAO, PetscTAOLineSearch*) 200 201 PetscErrorCode TaoALMMGetSubsolver(PetscTAO, PetscTAO*) 202 PetscErrorCode TaoALMMSetSubsolver(PetscTAO, PetscTAO) 203 PetscErrorCode TaoALMMGetType(PetscTAO, PetscTAOALMMType*) 204 PetscErrorCode TaoALMMSetType(PetscTAO, PetscTAOALMMType) 205 206 PetscErrorCode TaoBRGNGetSubsolver(PetscTAO, PetscTAO*) 207 PetscErrorCode TaoBRGNSetRegularizerObjectiveAndGradientRoutine(PetscTAO, PetscTaoRegularizerObjGrad*, void*) 208 PetscErrorCode TaoBRGNSetRegularizerHessianRoutine(PetscTAO, PetscMat, PetscTaoRegularizerHessian*, void*) 209 PetscErrorCode TaoBRGNSetRegularizerWeight(PetscTAO, PetscReal) 210 PetscErrorCode TaoBRGNSetL1SmoothEpsilon(PetscTAO, PetscReal) 211 PetscErrorCode TaoBRGNSetDictionaryMatrix(PetscTAO, PetscMat) 212 PetscErrorCode TaoBRGNGetDampingVector(PetscTAO, PetscVec*) 213 214 PetscErrorCode TaoPythonSetType(PetscTAO, char[]) 215 PetscErrorCode TaoPythonGetType(PetscTAO, char*[]) 216 217 ctypedef const char* PetscTAOLineSearchType "TaoLineSearchType" 218 PetscTAOLineSearchType TAOLINESEARCHUNIT 219 PetscTAOLineSearchType TAOLINESEARCHARMIJO 220 PetscTAOLineSearchType TAOLINESEARCHOWARMIJO 221 PetscTAOLineSearchType TAOLINESEARCHGPCG 222 PetscTAOLineSearchType TAOLINESEARCHMT 223 PetscTAOLineSearchType TAOLINESEARCHIPM 224 225 ctypedef enum PetscTAOLineSearchConvergedReason "TaoLineSearchConvergedReason": 226 # failed 227 TAOLINESEARCH_FAILED_INFORNAN 228 TAOLINESEARCH_FAILED_BADPARAMETER 229 TAOLINESEARCH_FAILED_ASCENT 230 # continue 231 TAOLINESEARCH_CONTINUE_ITERATING 232 # success 233 TAOLINESEARCH_SUCCESS 234 TAOLINESEARCH_SUCCESS_USER 235 # halted 236 TAOLINESEARCH_HALTED_OTHER 237 TAOLINESEARCH_HALTED_MAXFCN 238 TAOLINESEARCH_HALTED_UPPERBOUND 239 TAOLINESEARCH_HALTED_LOWERBOUND 240 TAOLINESEARCH_HALTED_RTOL 241 TAOLINESEARCH_HALTED_USER 242 243 ctypedef PetscErrorCode PetscTaoLineSearchObjective(PetscTAOLineSearch, PetscVec, PetscReal*, void*) except PETSC_ERR_PYTHON 244 ctypedef PetscErrorCode PetscTaoLineSearchGradient(PetscTAOLineSearch, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 245 ctypedef PetscErrorCode PetscTaoLineSearchObjGrad(PetscTAOLineSearch, PetscVec, PetscReal*, PetscVec, void*) except PETSC_ERR_PYTHON 246 ctypedef PetscErrorCode PetscTaoLineSearchObjGTS(PetscTaoLineSearch, PetscVec, PetscVec, PetscReal*, PetscReal*, void*) except PETSC_ERR_PYTHON 247 248 PetscErrorCode TaoLineSearchCreate(MPI_Comm, PetscTAOLineSearch*) 249 PetscErrorCode TaoLineSearchDestroy(PetscTAOLineSearch*) 250 PetscErrorCode TaoLineSearchView(PetscTAOLineSearch, PetscViewer) 251 PetscErrorCode TaoLineSearchSetType(PetscTAOLineSearch, PetscTAOLineSearchType) 252 PetscErrorCode TaoLineSearchGetType(PetscTAOLineSearch, PetscTAOLineSearchType*) 253 PetscErrorCode TaoLineSearchSetOptionsPrefix(PetscTAOLineSearch, char[]) 254 PetscErrorCode TaoLineSearchGetOptionsPrefix(PetscTAOLineSearch, char*[]) 255 PetscErrorCode TaoLineSearchSetFromOptions(PetscTAOLineSearch) 256 PetscErrorCode TaoLineSearchSetUp(PetscTAOLineSearch) 257 PetscErrorCode TaoLineSearchUseTaoRoutines(PetscTAOLineSearch, PetscTAO) 258 PetscErrorCode TaoLineSearchSetObjectiveRoutine(PetscTAOLineSearch, PetscTaoLineSearchObjective, void*) 259 PetscErrorCode TaoLineSearchSetGradientRoutine(PetscTAOLineSearch, PetscTaoLineSearchGradient, void*) 260 PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(PetscTAOLineSearch, PetscTaoLineSearchObjGrad, void*) 261 PetscErrorCode TaoLineSearchApply(PetscTAOLineSearch, PetscVec, PetscReal*, PetscVec, PetscVec, PetscReal*, PetscTAOLineSearchConvergedReason*) 262 PetscErrorCode TaoLineSearchSetInitialStepLength(PetscTAOLineSearch, PetscReal) 263 264# -------------------------------------------------------------------- 265 266cdef inline TAO ref_TAO(PetscTAO tao): 267 cdef TAO ob = <TAO> TAO() 268 ob.tao = tao 269 CHKERR(PetscINCREF(ob.obj)) 270 return ob 271 272# -------------------------------------------------------------------- 273 274cdef PetscErrorCode TAO_Objective(PetscTAO _tao, 275 PetscVec _x, PetscReal *_f, 276 void *ctx) except PETSC_ERR_PYTHON with gil: 277 278 cdef TAO tao = ref_TAO(_tao) 279 cdef Vec x = ref_Vec(_x) 280 context = tao.get_attr("__objective__") 281 if context is None and ctx != NULL: context = <object>ctx 282 assert context is not None and type(context) is tuple # sanity check 283 (objective, args, kargs) = context 284 retv = objective(tao, x, *args, **kargs) 285 _f[0] = asReal(retv) 286 return PETSC_SUCCESS 287 288cdef PetscErrorCode TAO_Residual(PetscTAO _tao, 289 PetscVec _x, PetscVec _r, 290 void *ctx) except PETSC_ERR_PYTHON with gil: 291 292 cdef TAO tao = ref_TAO(_tao) 293 cdef Vec x = ref_Vec(_x) 294 cdef Vec r = ref_Vec(_r) 295 context = tao.get_attr("__residual__") 296 if context is None and ctx != NULL: context = <object>ctx 297 assert context is not None and type(context) is tuple # sanity check 298 (residual, args, kargs) = context 299 residual(tao, x, r, *args, **kargs) 300 return PETSC_SUCCESS 301 302cdef PetscErrorCode TAO_Gradient(PetscTAO _tao, 303 PetscVec _x, PetscVec _g, 304 void *ctx) except PETSC_ERR_PYTHON with gil: 305 306 cdef TAO tao = ref_TAO(_tao) 307 cdef Vec x = ref_Vec(_x) 308 cdef Vec g = ref_Vec(_g) 309 context = tao.get_attr("__gradient__") 310 if context is None and ctx != NULL: context = <object>ctx 311 assert context is not None and type(context) is tuple # sanity check 312 (gradient, args, kargs) = context 313 gradient(tao, x, g, *args, **kargs) 314 return PETSC_SUCCESS 315 316cdef PetscErrorCode TAO_ObjGrad(PetscTAO _tao, 317 PetscVec _x, PetscReal *_f, PetscVec _g, 318 void *ctx) except PETSC_ERR_PYTHON with gil: 319 320 cdef TAO tao = ref_TAO(_tao) 321 cdef Vec x = ref_Vec(_x) 322 cdef Vec g = ref_Vec(_g) 323 context = tao.get_attr("__objgrad__") 324 if context is None and ctx != NULL: context = <object>ctx 325 assert context is not None and type(context) is tuple # sanity check 326 (objgrad, args, kargs) = context 327 retv = objgrad(tao, x, g, *args, **kargs) 328 _f[0] = asReal(retv) 329 return PETSC_SUCCESS 330 331 332cdef PetscErrorCode TAO_BRGNRegObjGrad(PetscTAO _tao, 333 PetscVec _x, PetscReal *_f, PetscVec _g, 334 void *ctx) except PETSC_ERR_PYTHON with gil: 335 336 cdef TAO tao = ref_TAO(_tao) 337 cdef Vec x = ref_Vec(_x) 338 cdef Vec g = ref_Vec(_g) 339 context = tao.get_attr("__brgnregobjgrad__") 340 if context is None and ctx != NULL: context = <object>ctx 341 assert context is not None and type(context) is tuple # sanity check 342 (objgrad, args, kargs) = context 343 retv = objgrad(tao, x, g, *args, **kargs) 344 _f[0] = asReal(retv) 345 return PETSC_SUCCESS 346 347cdef PetscErrorCode TAO_Constraints(PetscTAO _tao, 348 PetscVec _x, PetscVec _r, 349 void *ctx) except PETSC_ERR_PYTHON with gil: 350 351 cdef TAO tao = ref_TAO(_tao) 352 cdef Vec x = ref_Vec(_x) 353 cdef Vec r = ref_Vec(_r) 354 context = tao.get_attr("__constraints__") 355 if context is None and ctx != NULL: context = <object>ctx 356 assert context is not None and type(context) is tuple # sanity check 357 (constraints, args, kargs) = context 358 constraints(tao, x, r, *args, **kargs) 359 return PETSC_SUCCESS 360 361cdef PetscErrorCode TAO_VarBounds(PetscTAO _tao, 362 PetscVec _xl, PetscVec _xu, 363 void *ctx) except PETSC_ERR_PYTHON with gil: 364 365 cdef TAO tao = ref_TAO(_tao) 366 cdef Vec xl = ref_Vec(_xl) 367 cdef Vec xu = ref_Vec(_xu) 368 context = tao.get_attr("__varbounds__") 369 if context is None and ctx != NULL: context = <object>ctx 370 assert context is not None and type(context) is tuple # sanity check 371 (varbounds, args, kargs) = context 372 varbounds(tao, xl, xu, *args, **kargs) 373 return PETSC_SUCCESS 374 375cdef PetscErrorCode TAO_Hessian(PetscTAO _tao, 376 PetscVec _x, 377 PetscMat _H, 378 PetscMat _P, 379 void* ctx) except PETSC_ERR_PYTHON with gil: 380 cdef TAO tao = ref_TAO(_tao) 381 cdef Vec x = ref_Vec(_x) 382 cdef Mat H = ref_Mat(_H) 383 cdef Mat P = ref_Mat(_P) 384 context = tao.get_attr("__hessian__") 385 if context is None and ctx != NULL: context = <object>ctx 386 assert context is not None and type(context) is tuple # sanity check 387 (hessian, args, kargs) = context 388 hessian(tao, x, H, P, *args, **kargs) 389 return PETSC_SUCCESS 390 391cdef PetscErrorCode TAO_BRGNRegHessian(PetscTAO _tao, 392 PetscVec _x, 393 PetscMat _H, 394 void* ctx) except PETSC_ERR_PYTHON with gil: 395 cdef TAO tao = ref_TAO(_tao) 396 cdef Vec x = ref_Vec(_x) 397 cdef Mat H = ref_Mat(_H) 398 context = tao.get_attr("__brgnreghessian__") 399 if context is None and ctx != NULL: context = <object>ctx 400 assert context is not None and type(context) is tuple # sanity check 401 (hessian, args, kargs) = context 402 hessian(tao, x, H, *args, **kargs) 403 return PETSC_SUCCESS 404 405cdef PetscErrorCode TAO_Jacobian(PetscTAO _tao, 406 PetscVec _x, 407 PetscMat _J, 408 PetscMat _P, 409 void* ctx) except PETSC_ERR_PYTHON with gil: 410 cdef TAO tao = ref_TAO(_tao) 411 cdef Vec x = ref_Vec(_x) 412 cdef Mat J = ref_Mat(_J) 413 cdef Mat P = ref_Mat(_P) 414 context = tao.get_attr("__jacobian__") 415 if context is None and ctx != NULL: context = <object>ctx 416 assert context is not None and type(context) is tuple # sanity check 417 (jacobian, args, kargs) = context 418 jacobian(tao, x, J, P, *args, **kargs) 419 return PETSC_SUCCESS 420 421cdef PetscErrorCode TAO_JacobianResidual(PetscTAO _tao, 422 PetscVec _x, 423 PetscMat _J, 424 PetscMat _P, 425 void* ctx) except PETSC_ERR_PYTHON with gil: 426 cdef TAO tao = ref_TAO(_tao) 427 cdef Vec x = ref_Vec(_x) 428 cdef Mat J = ref_Mat(_J) 429 cdef Mat P = ref_Mat(_P) 430 context = tao.get_attr("__jacobian_residual__") 431 if context is None and ctx != NULL: context = <object>ctx 432 assert context is not None and type(context) is tuple # sanity check 433 (jacobian, args, kargs) = context 434 jacobian(tao, x, J, P, *args, **kargs) 435 return PETSC_SUCCESS 436 437cdef PetscErrorCode TAO_JacobianState(PetscTAO _tao, 438 PetscVec _x, 439 PetscMat _J, 440 PetscMat _Jp, 441 PetscMat _Ji, 442 void* ctx) except PETSC_ERR_PYTHON with gil: 443 cdef TAO tao = ref_TAO(_tao) 444 cdef Vec x = ref_Vec(_x) 445 cdef Mat J = ref_Mat(_J) 446 cdef Mat Jp = ref_Mat(_Jp) 447 cdef Mat Ji = ref_Mat(_Ji) 448 context = tao.get_attr("__jacobian_state__") 449 if context is None and ctx != NULL: context = <object>ctx 450 assert context is not None and type(context) is tuple # sanity check 451 (jacobian, args, kargs) = context 452 jacobian(tao, x, J, Jp, Ji, *args, **kargs) 453 return PETSC_SUCCESS 454 455cdef PetscErrorCode TAO_JacobianDesign(PetscTAO _tao, 456 PetscVec _x, 457 PetscMat _J, 458 void* ctx) except PETSC_ERR_PYTHON with gil: 459 cdef TAO tao = ref_TAO(_tao) 460 cdef Vec x = ref_Vec(_x) 461 cdef Mat J = ref_Mat(_J) 462 context = tao.get_attr("__jacobian_design__") 463 if context is None and ctx != NULL: context = <object>ctx 464 assert context is not None and type(context) is tuple # sanity check 465 (jacobian, args, kargs) = context 466 jacobian(tao, x, J, *args, **kargs) 467 return PETSC_SUCCESS 468 469cdef PetscErrorCode TAO_EqualityConstraints(PetscTAO _tao, 470 PetscVec _x, 471 PetscVec _c, 472 void* ctx) except PETSC_ERR_PYTHON with gil: 473 cdef TAO tao = ref_TAO(_tao) 474 cdef Vec x = ref_Vec(_x) 475 cdef Vec c = ref_Vec(_c) 476 context = tao.get_attr("__equality_constraints__") 477 if context is None and ctx != NULL: context = <object>ctx 478 assert context is not None and type(context) is tuple # sanity check 479 (f, args, kargs) = context 480 f(tao, x, c, *args, **kargs) 481 return PETSC_SUCCESS 482 483cdef PetscErrorCode TAO_InequalityConstraints(PetscTAO _tao, 484 PetscVec _x, 485 PetscVec _c, 486 void* ctx) except PETSC_ERR_PYTHON with gil: 487 cdef TAO tao = ref_TAO(_tao) 488 cdef Vec x = ref_Vec(_x) 489 cdef Vec c = ref_Vec(_c) 490 context = tao.get_attr("__inequality_constraints__") 491 if context is None and ctx != NULL: context = <object>ctx 492 assert context is not None and type(context) is tuple # sanity check 493 (f, args, kargs) = context 494 f(tao, x, c, *args, **kargs) 495 return PETSC_SUCCESS 496 497cdef PetscErrorCode TAO_JacobianEquality(PetscTAO _tao, 498 PetscVec _x, 499 PetscMat _J, 500 PetscMat _P, 501 void* ctx) except PETSC_ERR_PYTHON with gil: 502 cdef TAO tao = ref_TAO(_tao) 503 cdef Vec x = ref_Vec(_x) 504 cdef Mat J = ref_Mat(_J) 505 cdef Mat P = ref_Mat(_P) 506 context = tao.get_attr("__jacobian_equality__") 507 if context is None and ctx != NULL: context = <object>ctx 508 assert context is not None and type(context) is tuple # sanity check 509 (jacobian, args, kargs) = context 510 jacobian(tao, x, J, P, *args, **kargs) 511 return PETSC_SUCCESS 512 513cdef PetscErrorCode TAO_JacobianInequality(PetscTAO _tao, 514 PetscVec _x, 515 PetscMat _J, 516 PetscMat _P, 517 void* ctx) except PETSC_ERR_PYTHON with gil: 518 cdef TAO tao = ref_TAO(_tao) 519 cdef Vec x = ref_Vec(_x) 520 cdef Mat J = ref_Mat(_J) 521 cdef Mat P = ref_Mat(_P) 522 context = tao.get_attr("__jacobian_inequality__") 523 if context is None and ctx != NULL: context = <object>ctx 524 assert context is not None and type(context) is tuple # sanity check 525 (jacobian, args, kargs) = context 526 jacobian(tao, x, J, P, *args, **kargs) 527 return PETSC_SUCCESS 528 529# ctx is unused 530cdef PetscErrorCode TAO_Update( 531 PetscTAO _tao, 532 PetscInt its, 533 void* ctx) except PETSC_ERR_PYTHON with gil: 534 cdef TAO tao = ref_TAO(_tao) 535 (update, args, kargs) = tao.get_attr('__update__') 536 update(tao, toInt(its), *args, **kargs) 537 return PETSC_SUCCESS 538 539cdef PetscErrorCode TAO_Converged(PetscTAO _tao, 540 void* ctx) except PETSC_ERR_PYTHON with gil: 541 # call first the default convergence test 542 CHKERR(TaoDefaultConvergenceTest(_tao, NULL)) 543 # call next the user-provided convergence test 544 cdef TAO tao = ref_TAO(_tao) 545 (converged, args, kargs) = tao.get_attr('__converged__') 546 reason = converged(tao, *args, **kargs) 547 if reason is None: return PETSC_SUCCESS 548 # handle value of convergence reason 549 cdef PetscTAOConvergedReason creason = TAO_CONTINUE_ITERATING 550 if reason is False or reason == -1: 551 creason = TAO_DIVERGED_USER 552 elif reason is True or reason == 1: 553 creason = TAO_CONVERGED_USER 554 else: 555 creason = reason 556 assert creason >= TAO_DIVERGED_USER 557 assert creason <= TAO_CONVERGED_USER 558 CHKERR(TaoSetConvergedReason(_tao, creason)) 559 return PETSC_SUCCESS 560 561cdef PetscErrorCode TAO_Monitor(PetscTAO _tao, 562 void* ctx) except PETSC_ERR_PYTHON with gil: 563 cdef TAO tao = ref_TAO(_tao) 564 cdef object monitorlist = tao.get_attr('__monitor__') 565 if monitorlist is None: return PETSC_SUCCESS 566 for (monitor, args, kargs) in monitorlist: 567 monitor(tao, *args, **kargs) 568 return PETSC_SUCCESS 569 570# -------------------------------------------------------------------- 571 572cdef inline TAOLineSearch ref_TAOLS(PetscTAOLineSearch taols): 573 cdef TAOLineSearch ob = <TAOLineSearch> TAOLineSearch() 574 ob.taols = taols 575 CHKERR(PetscINCREF(ob.obj)) 576 return ob 577 578# -------------------------------------------------------------------- 579 580cdef PetscErrorCode TAOLS_Objective(PetscTAOLineSearch _ls, 581 PetscVec _x, 582 PetscReal *_f, 583 void *ctx) except PETSC_ERR_PYTHON with gil: 584 585 cdef TAOLineSearch ls = ref_TAOLS(_ls) 586 cdef Vec x = ref_Vec(_x) 587 (objective, args, kargs) = ls.get_attr("__objective__") 588 retv = objective(ls, x, *args, **kargs) 589 _f[0] = asReal(retv) 590 return PETSC_SUCCESS 591 592cdef PetscErrorCode TAOLS_Gradient(PetscTAOLineSearch _ls, 593 PetscVec _x, 594 PetscVec _g, 595 void *ctx) except PETSC_ERR_PYTHON with gil: 596 597 cdef TAOLineSearch ls = ref_TAOLS(_ls) 598 cdef Vec x = ref_Vec(_x) 599 cdef Vec g = ref_Vec(_g) 600 (gradient, args, kargs) = ls.get_attr("__gradient__") 601 gradient(ls, x, g, *args, **kargs) 602 return PETSC_SUCCESS 603 604 605cdef PetscErrorCode TAOLS_ObjGrad(PetscTAOLineSearch _ls, 606 PetscVec _x, 607 PetscReal *_f, 608 PetscVec _g, 609 void *ctx) except PETSC_ERR_PYTHON with gil: 610 611 cdef TAOLineSearch ls = ref_TAOLS(_ls) 612 cdef Vec x = ref_Vec(_x) 613 cdef Vec g = ref_Vec(_g) 614 (objgrad, args, kargs) = ls.get_attr("__objgrad__") 615 retv = objgrad(ls, x, g, *args, **kargs) 616 _f[0] = asReal(retv) 617 return PETSC_SUCCESS 618