1# cython: cdivision=True 2# cython: binding=False 3# cython: auto_pickle=False 4# cython: autotestdict=False 5# cython: warn.multiple_declarators=False 6 7# -------------------------------------------------------------------- 8 9cdef extern from "Python.h": 10 ctypedef struct PyObject 11 void Py_INCREF(PyObject*) 12 void Py_DECREF(PyObject*) 13 void Py_CLEAR(PyObject*) 14 object PyModule_New(const char[]) 15 bint PyModule_Check(object) 16 object PyImport_Import(object) 17 18# -------------------------------------------------------------------- 19 20cdef extern from * nogil: 21 PetscErrorCode PetscOptionsString(char[], char[], char[], char[], char[], size_t, PetscBool*) 22 23cdef extern from * nogil: # custom.h 24 PetscErrorCode PetscObjectComposedDataRegisterPy(PetscInt*) 25 PetscErrorCode PetscObjectComposedDataGetIntPy(PetscObject, PetscInt, PetscInt*, PetscBool*) 26 PetscErrorCode PetscObjectComposedDataSetIntPy(PetscObject, PetscInt, PetscInt) 27 28# -------------------------------------------------------------------- 29 30cdef char *FUNCT = NULL 31cdef char *fstack[1024] 32cdef int istack = 0 33 34cdef inline void FunctionBegin(char name[]) noexcept nogil: 35 global istack, fstack, FUNCT 36 FUNCT = name 37 fstack[istack] = FUNCT 38 istack += 1 39 if istack >= 1024: 40 istack = 0 41 return 42 43cdef inline PetscErrorCode FunctionEnd() noexcept nogil: 44 global istack, fstack, FUNCT 45 FUNCT = NULL 46 istack -= 1 47 if istack < 0: 48 istack = 1024 49 FUNCT = fstack[istack] 50 return PETSC_SUCCESS 51 52cdef PetscErrorCode PetscSETERR(PetscErrorCode ierr, const char msg[]) noexcept nogil: 53 global istack, fstack 54 istack = 0 55 fstack[istack] = NULL 56 return PetscERROR(PETSC_COMM_SELF, FUNCT, ierr, 57 PETSC_ERROR_INITIAL, msg, NULL) 58 59cdef PetscErrorCode UNSUPPORTED(const char msg[]) noexcept nogil: 60 return PetscERROR(PETSC_COMM_SELF, FUNCT, PETSC_ERR_USER, 61 PETSC_ERROR_INITIAL, b"method %s()", msg) 62 63# -------------------------------------------------------------------- 64 65cdef inline PetscInt getRef(void *pobj) noexcept nogil: 66 cdef PetscObject obj = <PetscObject>pobj 67 if obj == NULL: return 0 68 else: return obj.refct 69 70cdef inline void addRef(void *pobj) noexcept nogil: 71 cdef PetscObject obj = <PetscObject>pobj 72 if obj != NULL: obj.refct += 1 73 74cdef inline void delRef(void *pobj) noexcept nogil: 75 cdef PetscObject obj = <PetscObject>pobj 76 if obj != NULL: obj.refct -= 1 77 78cdef inline PetscObject newRef(void *pobj) noexcept nogil: 79 cdef PetscObject obj = <PetscObject>pobj 80 cdef int ierr = 0 81 if obj != NULL: 82 ierr = PetscObjectReference(obj) 83 if ierr: return NULL # XXX warning! 84 return obj 85 86cdef inline const char* getPrefix(void *pobj) noexcept nogil: 87 cdef PetscObject obj = <PetscObject>pobj 88 if obj == NULL: return NULL 89 return obj.prefix 90 91cdef inline int getCommSize(void *pobj) noexcept nogil: 92 cdef PetscObject obj = <PetscObject>pobj 93 if obj == NULL: return 0 94 cdef int size = 0 95 MPI_Comm_size(obj.comm, &size) 96 return size 97 98cdef inline Viewer Viewer_(PetscViewer p): 99 cdef Viewer ob = Viewer.__new__(Viewer) 100 ob.obj[0] = newRef(p) 101 return ob 102 103cdef inline IS IS_(PetscIS p): 104 cdef IS ob = IS.__new__(IS) 105 ob.obj[0] = newRef(p) 106 return ob 107 108cdef inline Vec Vec_(PetscVec p): 109 cdef Vec ob = Vec.__new__(Vec) 110 ob.obj[0] = newRef(p) 111 return ob 112 113cdef inline Mat Mat_(PetscMat p): 114 cdef Mat ob = Mat.__new__(Mat) 115 ob.obj[0] = newRef(p) 116 return ob 117 118cdef inline PC PC_(PetscPC p): 119 cdef PC ob = PC.__new__(PC) 120 ob.obj[0] = newRef(p) 121 return ob 122 123cdef inline KSP KSP_(PetscKSP p): 124 cdef KSP ob = KSP.__new__(KSP) 125 ob.obj[0] = newRef(p) 126 return ob 127 128cdef inline SNES SNES_(PetscSNES p): 129 cdef SNES ob = SNES.__new__(SNES) 130 ob.obj[0] = newRef(p) 131 return ob 132 133cdef inline TS TS_(PetscTS p): 134 cdef TS ob = TS.__new__(TS) 135 ob.obj[0] = newRef(p) 136 return ob 137 138cdef inline TAO TAO_(PetscTAO p): 139 cdef TAO ob = TAO.__new__(TAO) 140 ob.obj[0] = newRef(p) 141 return ob 142 143# -------------------------------------------------------------------- 144 145cdef object parse_url(object url): 146 path, name = url.rsplit(":", 1) 147 return (path, name) 148 149cdef dict module_cache = {} 150 151cdef object load_module(object path): 152 if path in module_cache: 153 return module_cache[path] 154 module = PyModule_New("__petsc__") 155 module.__file__ = path 156 module.__package__ = None 157 module_cache[path] = module 158 cdef object code = None 159 try: 160 with open(path, 'r') as source: 161 code = compile(source.read(), path, 'exec') 162 exec(code, module.__dict__) 163 except Exception: 164 del module_cache[path] 165 raise 166 return module 167 168# ----------------------------------------------------------------------------- 169 170 171@cython.internal 172cdef class _PyObj: 173 174 cdef object self 175 cdef bytes name 176 177 def __getattr__(self, attr): 178 return getattr(self.self, attr, None) 179 180 cdef int setcontext(self, void *ctx, Object base) except -1: 181 # 182 if ctx == <void*>self.self: 183 return 0 184 # 185 cdef object destroy = self.destroy 186 if destroy is not None: 187 destroy(base) 188 destroy = None 189 # 190 if ctx == NULL: 191 self.self = None 192 self.name = None 193 return 0 194 # 195 self.self = <object>ctx 196 self.name = None 197 cdef object create = self.create 198 if create is not None: 199 create(base) 200 create = None 201 return 0 202 203 cdef int getcontext(self, void **ctx) except -1: 204 if ctx == NULL: return 0 205 if self.self is not None: 206 ctx[0] = <void*> self.self 207 else: 208 ctx[0] = NULL 209 return 0 210 211 cdef int setname(self, const char name[]) except -1: 212 if name != NULL and name[0] != 0: 213 self.name = name 214 else: 215 self.name = None 216 return 0 217 218 cdef char* getname(self) except? NULL: 219 if self.self is None: 220 return NULL 221 if self.name is not None: 222 return self.name 223 cdef ctx = self.self 224 cdef name = None 225 if PyModule_Check(ctx): 226 name = getattr(ctx, '__name__', None) 227 else: 228 modname = getattr(ctx, '__module__', None) 229 clsname = None 230 cls = getattr(ctx, '__class__', None) 231 if cls: 232 clsname = getattr(cls, '__name__', None) 233 if not modname: 234 modname = getattr(cls, '__module__', None) 235 if modname and clsname: 236 name = modname + '.' + clsname 237 elif clsname: 238 name = clsname 239 elif modname: 240 name = modname 241 if name is not None: 242 self.name = name.encode() 243 if self.name is not None: 244 return self.name 245 return NULL 246 247cdef createcontext(const char name_p[]): 248 if name_p == NULL: return None 249 cdef name = bytes2str(name_p) 250 cdef mod 251 cdef path 252 cdef modname=None 253 cdef cls 254 cdef attr 255 cdef clsname=None 256 # path/to/filename.py:{function|class} 257 if ':' in name: 258 path, attr = parse_url(name) 259 mod = load_module(path) 260 if attr: 261 cls = getattr(mod, attr) 262 return cls() 263 else: 264 return mod 265 # package.module[.{function|class}] 266 if '.' in name: 267 modname, clsname = name.rsplit('.', 1) 268 mod = PyImport_Import(modname) 269 if hasattr(mod, clsname): 270 cls = getattr(mod, clsname) 271 if not PyModule_Check(cls): 272 return cls() 273 # package[.module] 274 mod = PyImport_Import(name) 275 return mod 276 277cdef int viewcontext(_PyObj ctx, PetscViewer viewer) except -1: 278 cdef PetscBool isascii = PETSC_FALSE, isstring = PETSC_FALSE 279 CHKERR(PetscObjectTypeCompare(<PetscObject>viewer, PETSCVIEWERASCII, &isascii)) 280 CHKERR(PetscObjectTypeCompare(<PetscObject>viewer, PETSCVIEWERSTRING, &isstring)) 281 cdef char *name = ctx.getname() 282 if isascii: 283 if name == NULL: name = b"unknown/not yet set" 284 CHKERR(PetscViewerASCIIPrintf(viewer, b" Python: %s\n", name)) 285 if isstring: 286 if name == NULL: name = b"<unknown>" 287 CHKERR(PetscViewerStringSPrintf(viewer, "%s", name)) 288 return 0 289 290# -------------------------------------------------------------------- 291 292cdef extern from * nogil: 293 struct _PetscViewerOps: 294 PetscErrorCode (*destroy)(PetscViewer) except PETSC_ERR_PYTHON 295 PetscErrorCode (*setup)(PetscViewer) except PETSC_ERR_PYTHON 296 PetscErrorCode (*flush)(PetscViewer) except PETSC_ERR_PYTHON 297 PetscErrorCode (*setfromoptions)(PetscViewer, PetscOptionItems) except PETSC_ERR_PYTHON 298 PetscErrorCode (*view)(PetscViewer, PetscViewer) except PETSC_ERR_PYTHON 299 ctypedef _PetscViewerOps *PetscViewerOps 300 struct _p_PetscViewer: 301 void *data 302 PetscViewerOps ops 303 304 305@cython.internal 306cdef class _PyVwr(_PyObj): 307 308 cdef bytes filename 309 310 cdef int setfilename(self, const char name[]) except -1: 311 if name != NULL and name[0] != 0: 312 self.filename = name 313 else: 314 self.filename = None 315 return 0 316 317 cdef char* getfilename(self) except? NULL: 318 if self.self is None: 319 return NULL 320 if self.filename is not None: 321 return self.filename 322 return NULL 323 324cdef inline _PyVwr PyVwr(PetscViewer viewer): 325 if viewer != NULL and viewer.data != NULL: 326 return <_PyVwr>viewer.data 327 else: 328 return _PyVwr.__new__(_PyVwr) 329 330cdef public PetscErrorCode PetscViewerPythonGetContext(PetscViewer viewer, void **ctx) \ 331 except PETSC_ERR_PYTHON: 332 FunctionBegin(b"PetscViewerPythonGetContext") 333 PyVwr(viewer).getcontext(ctx) 334 return FunctionEnd() 335 336cdef public PetscErrorCode PetscViewerPythonSetContext(PetscViewer viewer, void *ctx) \ 337 except PETSC_ERR_PYTHON: 338 FunctionBegin(b"PetscViewerPythonSetContext") 339 PyVwr(viewer).setcontext(ctx, Viewer_(viewer)) 340 return FunctionEnd() 341 342cdef PetscErrorCode PetscViewerPythonSetType_PYTHON(PetscViewer viewer, const char *name) \ 343 except PETSC_ERR_PYTHON with gil: 344 FunctionBegin(b"PetscViewerPythonSetType_PYTHON") 345 if name == NULL: return FunctionEnd() # XXX 346 cdef object ctx = createcontext(name) 347 PetscViewerPythonSetContext(viewer, <void*>ctx) 348 PyVwr(viewer).setname(name) 349 return FunctionEnd() 350 351cdef PetscErrorCode PetscViewerPythonGetType_PYTHON(PetscViewer viewer, const char *name[]) \ 352 except PETSC_ERR_PYTHON with gil: 353 FunctionBegin(b"PetscViewerPythonGetType_PYTHON") 354 name[0] = PyVwr(viewer).getname() 355 return FunctionEnd() 356 357cdef PetscErrorCode PetscViewerPythonSetFilename_PYTHON(PetscViewer viewer, const char *name) \ 358 except PETSC_ERR_PYTHON with gil: 359 FunctionBegin(b"PetscViewerPythonSetFilename_PYTHON") 360 PyVwr(viewer).setfilename(name) 361 return FunctionEnd() 362 363cdef PetscErrorCode PetscViewerPythonGetFilename_PYTHON(PetscViewer viewer, const char *name[]) \ 364 except PETSC_ERR_PYTHON with gil: 365 FunctionBegin(b"PetscViewerPythonGetFilename_PYTHON") 366 name[0] = PyVwr(viewer).getfilename() 367 return FunctionEnd() 368 369cdef PetscErrorCode PetscViewerPythonViewObject_PYTHON(PetscViewer viewer, PetscObject obj) \ 370 except PETSC_ERR_PYTHON with gil: 371 FunctionBegin(b"PetscViewerPythonViewObject_PYTHON") 372 cdef viewObject = PyVwr(viewer).viewObject 373 cdef Object pobj = PyPetscObject_New(obj) 374 if viewObject is not None: 375 viewObject(Viewer_(viewer), pobj) 376 return FunctionEnd() 377 378cdef PetscErrorCode PetscViewerCreate_Python( 379 PetscViewer viewer, 380 ) except PETSC_ERR_PYTHON with gil: 381 FunctionBegin(b"PetscViewerCreate_Python") 382 cdef PetscViewerOps ops = viewer.ops 383 ops.destroy = PetscViewerDestroy_Python 384 ops.view = PetscViewerView_Python 385 ops.setup = PetscViewerSetUp_Python 386 ops.setfromoptions = PetscViewerSetFromOptions_Python 387 ops.flush = PetscViewerFlush_Python 388 # 389 CHKERR(PetscObjectComposeFunction( 390 <PetscObject>viewer, b"PetscViewerPythonSetType_C", 391 <PetscVoidFunction>PetscViewerPythonSetType_PYTHON)) 392 CHKERR(PetscObjectComposeFunction( 393 <PetscObject>viewer, b"PetscViewerPythonGetType_C", 394 <PetscVoidFunction>PetscViewerPythonGetType_PYTHON)) 395 CHKERR(PetscObjectComposeFunction( 396 <PetscObject>viewer, b"PetscViewerFileSetName_C", 397 <PetscVoidFunction>PetscViewerPythonSetFilename_PYTHON)) 398 CHKERR(PetscObjectComposeFunction( 399 <PetscObject>viewer, b"PetscViewerFileGetName_C", 400 <PetscVoidFunction>PetscViewerPythonGetFilename_PYTHON)) 401 CHKERR(PetscObjectComposeFunction( 402 <PetscObject>viewer, b"PetscViewerPythonViewObject_C", 403 <PetscVoidFunction>PetscViewerPythonViewObject_PYTHON)) 404 # 405 cdef ctx = PyVwr(NULL) 406 viewer.data = <void*> ctx 407 Py_INCREF(<PyObject*>viewer.data) 408 return FunctionEnd() 409 410cdef inline PetscErrorCode PetscViewerDestroy_Python_inner( 411 PetscViewer viewer, 412 ) except PETSC_ERR_PYTHON with gil: 413 try: 414 addRef(viewer) 415 PetscViewerPythonSetContext(viewer, NULL) 416 finally: 417 delRef(viewer) 418 Py_DECREF(<PyObject*>viewer.data) 419 viewer.data = NULL 420 return PETSC_SUCCESS 421 422cdef PetscErrorCode PetscViewerDestroy_Python( 423 PetscViewer viewer, 424 ) except PETSC_ERR_PYTHON nogil: 425 FunctionBegin(b"PetscViewerDestroy_Python") 426 CHKERR(PetscObjectComposeFunction( 427 <PetscObject>viewer, b"PetscViewerPythonSetType_C", 428 <PetscVoidFunction>NULL)) 429 CHKERR(PetscObjectComposeFunction( 430 <PetscObject>viewer, b"PetscViewerPythonGetType_C", 431 <PetscVoidFunction>NULL)) 432 CHKERR(PetscObjectComposeFunction( 433 <PetscObject>viewer, b"PetscViewerFileSetName_C", 434 <PetscVoidFunction>NULL)) 435 CHKERR(PetscObjectComposeFunction( 436 <PetscObject>viewer, b"PetscViewerFileGetName_C", 437 <PetscVoidFunction>NULL)) 438 CHKERR(PetscObjectComposeFunction( 439 <PetscObject>viewer, b"PetscViewerPythonViewObject_C", 440 <PetscVoidFunction>NULL)) 441 # 442 if Py_IsInitialized(): PetscViewerDestroy_Python_inner(viewer) 443 return FunctionEnd() 444 445cdef PetscErrorCode PetscViewerSetUp_Python( 446 PetscViewer viewer, 447 ) except PETSC_ERR_PYTHON with gil: 448 FunctionBegin(b"PetscViewerSetUp_Python") 449 cdef char name[2048] 450 cdef PetscBool found = PETSC_FALSE 451 if PyVwr(viewer).self is None: 452 CHKERR(PetscOptionsGetString(NULL, 453 getPrefix(viewer), b"-viewer_python_type", 454 name, sizeof(name), &found)) 455 if found and name[0]: 456 CHKERR(PetscViewerPythonSetType_PYTHON(viewer, name)) 457 cdef char fname[2048] 458 cdef const char *fdefault = NULL 459 CHKERR(PetscViewerFileGetName(viewer, &fdefault)) 460 if not fdefault: 461 CHKERR(PetscOptionsGetString(NULL, 462 getPrefix(viewer), b"-viewer_python_filename", 463 fname, sizeof(fname), &found)) 464 if found and fname[0]: 465 CHKERR(PetscViewerPythonSetFilename_PYTHON(viewer, fname)) 466 if PyVwr(viewer).self is None: 467 return PetscSETERR(PETSC_ERR_USER, 468 "Python context not set, call one of \n" 469 " * PetscViewerPythonSetType(viewer, \"[package.]module.class\")\n" 470 " * PetscViewerSetFromOptions(viewer) and pass option " 471 "-viewer_python_type [package.]module.class") 472 # 473 cdef setUp = PyVwr(viewer).setUp 474 if setUp is not None: 475 setUp(Viewer_(viewer)) 476 return FunctionEnd() 477 478cdef PetscErrorCode PetscViewerSetFromOptions_Python( 479 PetscViewer viewer, 480 PetscOptionItems PetscOptionsObject, 481 ) except PETSC_ERR_PYTHON with gil: 482 FunctionBegin(b"PetscViewerSetFromOptions_Python") 483 cdef char name[2048], *defval = PyVwr(viewer).getname() 484 cdef PetscBool found = PETSC_FALSE 485 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 486 CHKERR(PetscOptionsString( 487 b"-viewer_python_type", b"Python [package.]module[.{class|function}]", 488 b"PetscViewerPythonSetType", defval, name, sizeof(name), &found)); <void>opts 489 if found and name[0]: 490 CHKERR(PetscViewerPythonSetType_PYTHON(viewer, name)) 491 # 492 cdef char fname[2048], *fdefval = PyVwr(viewer).getfilename() 493 CHKERR(PetscOptionsString( 494 b"-viewer_python_filename", b"Output filename for viewer", 495 b"PetscViewerPythonSetFilename", fdefval, fname, sizeof(fname), &found)); <void>opts 496 if found and fname[0]: 497 CHKERR(PetscViewerPythonSetFilename_PYTHON(viewer, fname)) 498 cdef setFromOptions = PyVwr(viewer).setFromOptions 499 if setFromOptions is not None: 500 setFromOptions(Viewer_(viewer)) 501 return FunctionEnd() 502 503cdef PetscErrorCode PetscViewerView_Python( 504 PetscViewer viewer, 505 PetscViewer vwr, 506 ) except PETSC_ERR_PYTHON with gil: 507 FunctionBegin(b"PetscViewerView_Python") 508 viewcontext(PyVwr(viewer), vwr) 509 cdef view = PyVwr(viewer).view 510 if view is not None: 511 view(Viewer_(viewer), Viewer_(vwr)) 512 return FunctionEnd() 513 514cdef PetscErrorCode PetscViewerFlush_Python( 515 PetscViewer viewer, 516 ) except PETSC_ERR_PYTHON with gil: 517 FunctionBegin(b"PetscViewerFlush_Python") 518 cdef flush = PyVwr(viewer).flush 519 if flush is not None: 520 flush(Viewer_(viewer)) 521 return FunctionEnd() 522# -------------------------------------------------------------------- 523 524cdef extern from * nogil: 525 struct _MatOps: 526 PetscErrorCode (*destroy)(PetscMat) except PETSC_ERR_PYTHON 527 PetscErrorCode (*setfromoptions)(PetscMat, PetscOptionItems) except PETSC_ERR_PYTHON 528 PetscErrorCode (*view)(PetscMat, PetscViewer) except PETSC_ERR_PYTHON 529 PetscErrorCode (*duplicate)(PetscMat, PetscMatDuplicateOption, PetscMat*) except PETSC_ERR_PYTHON 530 PetscErrorCode (*copy)(PetscMat, PetscMat, PetscMatStructure) except PETSC_ERR_PYTHON 531 PetscErrorCode (*createsubmatrix)(PetscMat, PetscIS, PetscIS, PetscMatReuse, PetscMat*) except PETSC_ERR_PYTHON 532 PetscErrorCode (*setoption)(PetscMat, PetscMatOption, PetscBool) except PETSC_ERR_PYTHON 533 PetscErrorCode (*setup)(PetscMat) except PETSC_ERR_PYTHON 534 PetscErrorCode (*assemblybegin)(PetscMat, PetscMatAssemblyType) except PETSC_ERR_PYTHON 535 PetscErrorCode (*assemblyend)(PetscMat, PetscMatAssemblyType) except PETSC_ERR_PYTHON 536 PetscErrorCode (*zeroentries)(PetscMat) except PETSC_ERR_PYTHON 537 PetscErrorCode (*zerorowscolumns)(PetscMat, PetscInt, PetscInt*, PetscScalar, PetscVec, PetscVec) except PETSC_ERR_PYTHON 538 PetscErrorCode (*scale)(PetscMat, PetscScalar) except PETSC_ERR_PYTHON 539 PetscErrorCode (*shift)(PetscMat, PetscScalar) except PETSC_ERR_PYTHON 540 PetscErrorCode (*sor)(PetscMat, PetscVec, PetscReal, PetscMatSORType, PetscReal, PetscInt, PetscInt, PetscVec) except PETSC_ERR_PYTHON 541 PetscErrorCode (*getvecs)(PetscMat, PetscVec*, PetscVec*) except PETSC_ERR_PYTHON 542 PetscErrorCode (*mult)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 543 PetscErrorCode (*multtranspose)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 544 PetscErrorCode (*multhermitian"multhermitiantranspose")(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 545 PetscErrorCode (*multadd)(PetscMat, PetscVec, PetscVec, PetscVec) except PETSC_ERR_PYTHON 546 PetscErrorCode (*multtransposeadd)(PetscMat, PetscVec, PetscVec, PetscVec) except PETSC_ERR_PYTHON 547 PetscErrorCode (*multhermitianadd"multhermitiantransposeadd")(PetscMat, PetscVec, PetscVec, PetscVec) except PETSC_ERR_PYTHON 548 PetscErrorCode (*multdiagonalblock)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 549 PetscErrorCode (*solve)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 550 PetscErrorCode (*solvetranspose)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 551 PetscErrorCode (*solveadd)(PetscMat, PetscVec, PetscVec, PetscVec) except PETSC_ERR_PYTHON 552 PetscErrorCode (*solvetransposeadd)(PetscMat, PetscVec, PetscVec, PetscVec) except PETSC_ERR_PYTHON 553 PetscErrorCode (*getdiagonal)(PetscMat, PetscVec) except PETSC_ERR_PYTHON 554 PetscErrorCode (*setdiagonal"diagonalset")(PetscMat, PetscVec, PetscInsertMode) except PETSC_ERR_PYTHON 555 PetscErrorCode (*diagonalscale)(PetscMat, PetscVec, PetscVec) except PETSC_ERR_PYTHON 556 PetscErrorCode (*norm)(PetscMat, PetscNormType, PetscReal*) except PETSC_ERR_PYTHON 557 PetscErrorCode (*realpart)(PetscMat) except PETSC_ERR_PYTHON 558 PetscErrorCode (*imagpart"imaginarypart")(PetscMat) except PETSC_ERR_PYTHON 559 PetscErrorCode (*conjugate)(PetscMat) except PETSC_ERR_PYTHON 560 PetscErrorCode (*getdiagonalblock)(PetscMat, PetscMat*) except PETSC_ERR_PYTHON 561 PetscErrorCode (*productsetfromoptions)(PetscMat) except PETSC_ERR_PYTHON 562 PetscErrorCode (*productsymbolic)(PetscMat) except PETSC_ERR_PYTHON 563 PetscErrorCode (*productnumeric)(PetscMat) except PETSC_ERR_PYTHON 564 PetscErrorCode (*hasoperation)(PetscMat, PetscMatOperation, PetscBool*) except PETSC_ERR_PYTHON 565 ctypedef _MatOps *MatOps 566 ctypedef struct Mat_Product: 567 void *data 568 struct _p_Mat: 569 void *data 570 MatOps ops 571 PetscBool assembled 572 PetscBool preallocated 573 PetscLayout rmap, cmap 574 Mat_Product *product 575 576 577@cython.internal 578cdef class _PyMat(_PyObj): pass 579cdef inline _PyMat PyMat(PetscMat mat): 580 if mat != NULL and mat.data != NULL: 581 return <_PyMat>mat.data 582 else: 583 return _PyMat.__new__(_PyMat) 584 585cdef public PetscErrorCode MatPythonGetContext(PetscMat mat, void **ctx) \ 586 except PETSC_ERR_PYTHON: 587 FunctionBegin(b"MatPythonGetContext") 588 PyMat(mat).getcontext(ctx) 589 return FunctionEnd() 590 591cdef public PetscErrorCode MatPythonSetContext(PetscMat mat, void *ctx) \ 592 except PETSC_ERR_PYTHON: 593 FunctionBegin(b"MatPythonSetContext") 594 PyMat(mat).setcontext(ctx, Mat_(mat)) 595 mat.preallocated = PETSC_FALSE 596 return FunctionEnd() 597 598cdef PetscErrorCode MatPythonSetType_PYTHON(PetscMat mat, const char *name) \ 599 except PETSC_ERR_PYTHON with gil: 600 FunctionBegin(b"MatPythonSetType_PYTHON") 601 if name == NULL: return FunctionEnd() # XXX 602 cdef object ctx = createcontext(name) 603 MatPythonSetContext(mat, <void*>ctx) 604 PyMat(mat).setname(name) 605 return FunctionEnd() 606 607cdef PetscErrorCode MatPythonGetType_PYTHON(PetscMat mat, const char *name[]) \ 608 except PETSC_ERR_PYTHON with gil: 609 FunctionBegin(b"MatPythonGetType_PYTHON") 610 name[0] = PyMat(mat).getname() 611 return FunctionEnd() 612 613# FIXME: view and setFromOptions? 614cdef dict dMatOps = { 615 3 : 'mult', 616 4 : 'multAdd', 617 5 : 'multTranspose', 618 6 : 'multTransposeAdd', 619 7 : 'solve', 620 8 : 'solveAdd', 621 9 : 'solveTranspose', 622 10 : 'solveTransposeAdd', 623 13 : 'SOR', 624 17 : 'getDiagonal', 625 18 : 'diagonalScale', 626 19 : 'norm', 627 23 : 'zeroEntries', 628 32 : 'getDiagonalBlock', 629 34 : 'duplicate', 630 43 : 'copy', 631 45 : 'scale', 632 46 : 'shift', 633 47 : 'setDiagonal', 634 48 : 'zeroRowsColumns', 635 59 : 'createSubMatrix', 636 83 : 'getVecs', # FIXME -> createVecs 637 93 : 'conjugate', 638 96 : 'realPart', 639 97 : 'imagPart', 640 109 : 'multDiagonalBlock', 641 111 : 'multHermitian', 642 112 : 'multHermitianAdd', 643 } 644 645cdef PetscErrorCode MatCreate_Python( 646 PetscMat mat, 647 ) except PETSC_ERR_PYTHON with gil: 648 FunctionBegin(b"MatCreate_Python") 649 cdef MatOps ops = mat.ops 650 ops.destroy = MatDestroy_Python 651 ops.setfromoptions = MatSetFromOptions_Python 652 ops.view = MatView_Python 653 ops.duplicate = MatDuplicate_Python 654 ops.copy = MatCopy_Python 655 ops.createsubmatrix = MatCreateSubMatrix_Python 656 ops.setoption = MatSetOption_Python 657 ops.setup = MatSetUp_Python 658 ops.assemblybegin = MatAssemblyBegin_Python 659 ops.assemblyend = MatAssemblyEnd_Python 660 ops.zeroentries = MatZeroEntries_Python 661 ops.zerorowscolumns = MatZeroRowsColumns_Python 662 ops.scale = MatScale_Python 663 ops.shift = MatShift_Python 664 ops.getvecs = MatCreateVecs_Python 665 ops.mult = MatMult_Python 666 ops.sor = MatSOR_Python 667 ops.multtranspose = MatMultTranspose_Python 668 ops.multhermitian = MatMultHermitian_Python 669 ops.multadd = MatMultAdd_Python 670 ops.multtransposeadd = MatMultTransposeAdd_Python 671 ops.multhermitianadd = MatMultHermitianAdd_Python 672 ops.multdiagonalblock = MatMultDiagonalBlock_Python 673 ops.solve = MatSolve_Python 674 ops.solvetranspose = MatSolveTranspose_Python 675 ops.solveadd = MatSolveAdd_Python 676 ops.solvetransposeadd = MatSolveTransposeAdd_Python 677 ops.getdiagonal = MatGetDiagonal_Python 678 ops.setdiagonal = MatSetDiagonal_Python 679 ops.diagonalscale = MatDiagonalScale_Python 680 ops.norm = MatNorm_Python 681 ops.realpart = MatRealPart_Python 682 ops.imagpart = MatImagPart_Python 683 ops.conjugate = MatConjugate_Python 684 ops.hasoperation = MatHasOperation_Python 685 ops.getdiagonalblock = MatGetDiagonalBlock_Python 686 687 ops.productsetfromoptions = MatProductSetFromOptions_Python 688 689 # 690 mat.assembled = PETSC_TRUE # XXX 691 mat.preallocated = PETSC_FALSE # XXX 692 # 693 CHKERR(PetscObjectComposeFunction( 694 <PetscObject>mat, b"MatPythonSetType_C", 695 <PetscVoidFunction>MatPythonSetType_PYTHON)) 696 CHKERR(PetscObjectComposeFunction( 697 <PetscObject>mat, b"MatPythonGetType_C", 698 <PetscVoidFunction>MatPythonGetType_PYTHON)) 699 CHKERR(PetscObjectComposeFunction( 700 <PetscObject>mat, b"MatProductSetFromOptions_anytype_C", 701 <PetscVoidFunction>MatProductSetFromOptions_Python)) 702 CHKERR(PetscObjectChangeTypeName( 703 <PetscObject>mat, MATPYTHON)) 704 # 705 cdef ctx = PyMat(NULL) 706 mat.data = <void*> ctx 707 Py_INCREF(<PyObject*>mat.data) 708 return FunctionEnd() 709 710cdef inline PetscErrorCode MatDestroy_Python_inner( 711 PetscMat mat, 712 ) except PETSC_ERR_PYTHON with gil: 713 try: 714 addRef(mat) 715 MatPythonSetContext(mat, NULL) 716 finally: 717 delRef(mat) 718 Py_DECREF(<PyObject*>mat.data) 719 mat.data = NULL 720 return PETSC_SUCCESS 721 722cdef PetscErrorCode MatDestroy_Python( 723 PetscMat mat, 724 ) except PETSC_ERR_PYTHON nogil: 725 FunctionBegin(b"MatDestroy_Python") 726 CHKERR(PetscObjectComposeFunction( 727 <PetscObject>mat, b"MatPythonSetType_C", 728 <PetscVoidFunction>NULL)) 729 CHKERR(PetscObjectComposeFunction( 730 <PetscObject>mat, b"MatPythonGetType_C", 731 <PetscVoidFunction>NULL)) 732 CHKERR(PetscObjectComposeFunction( 733 <PetscObject>mat, b"MatProductSetFromOptions_anytype_C", 734 <PetscVoidFunction>NULL)) 735 CHKERR(PetscObjectChangeTypeName( 736 <PetscObject>mat, NULL)) 737 738 if Py_IsInitialized(): MatDestroy_Python_inner(mat) 739 return FunctionEnd() 740 741cdef PetscErrorCode MatSetFromOptions_Python( 742 PetscMat mat, 743 PetscOptionItems PetscOptionsObject, 744 ) except PETSC_ERR_PYTHON with gil: 745 FunctionBegin(b"MatSetFromOptions_Python") 746 cdef char name[2048], *defval = PyMat(mat).getname() 747 cdef PetscBool found = PETSC_FALSE 748 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 749 CHKERR(PetscOptionsString( 750 b"-mat_python_type", b"Python [package.]module[.{class|function}]", 751 b"MatPythonSetType", defval, name, sizeof(name), &found)); <void>opts 752 if found and name[0]: 753 CHKERR(MatPythonSetType_PYTHON(mat, name)) 754 # 755 cdef setFromOptions = PyMat(mat).setFromOptions 756 if setFromOptions is not None: 757 setFromOptions(Mat_(mat)) 758 return FunctionEnd() 759 760cdef PetscErrorCode MatView_Python( 761 PetscMat mat, 762 PetscViewer vwr, 763 ) except PETSC_ERR_PYTHON with gil: 764 FunctionBegin(b"MatView_Python") 765 viewcontext(PyMat(mat), vwr) 766 cdef view = PyMat(mat).view 767 if view is not None: 768 view(Mat_(mat), Viewer_(vwr)) 769 return FunctionEnd() 770 771cdef PetscErrorCode MatDuplicate_Python( 772 PetscMat mat, 773 PetscMatDuplicateOption op, 774 PetscMat* out, 775 ) except PETSC_ERR_PYTHON with gil: 776 FunctionBegin(b"MatDuplicate_Python") 777 cdef duplicate = PyMat(mat).duplicate 778 if duplicate is None: return UNSUPPORTED(b"duplicate") 779 cdef Mat m = duplicate(Mat_(mat), <long>op) 780 out[0] = m.mat; m.mat = NULL 781 return FunctionEnd() 782 783cdef PetscErrorCode MatCopy_Python( 784 PetscMat mat, 785 PetscMat out, 786 PetscMatStructure op, 787 ) except PETSC_ERR_PYTHON with gil: 788 FunctionBegin(b"MatCopy_Python") 789 cdef copy = PyMat(mat).copy 790 if copy is None: return UNSUPPORTED(b"copy") 791 copy(Mat_(mat), Mat_(out), <long>op) 792 return FunctionEnd() 793 794cdef PetscErrorCode MatGetDiagonalBlock_Python( 795 PetscMat mat, 796 PetscMat *out 797 ) except PETSC_ERR_PYTHON with gil: 798 FunctionBegin(b"MatGetDiagonalBlock_Python") 799 cdef getDiagonalBlock = PyMat(mat).getDiagonalBlock 800 if getDiagonalBlock is None: 801 try: 802 mat.ops.getdiagonalblock = NULL 803 CHKERR(MatGetDiagonalBlock(mat, out)) 804 finally: 805 mat.ops.getdiagonalblock = MatGetDiagonalBlock_Python 806 return FunctionEnd() 807 cdef Mat sub = getDiagonalBlock(Mat_(mat)) 808 if sub is not None: out[0] = sub.mat 809 return FunctionEnd() 810 811cdef PetscErrorCode MatCreateSubMatrix_Python( 812 PetscMat mat, 813 PetscIS row, 814 PetscIS col, 815 PetscMatReuse op, 816 PetscMat *out, 817 ) except PETSC_ERR_PYTHON with gil: 818 FunctionBegin(b"MatCreateSubMatrix_Python") 819 cdef createSubMatrix = PyMat(mat).createSubMatrix 820 if createSubMatrix is None: 821 try: 822 mat.ops.createsubmatrix = NULL 823 CHKERR(MatCreateSubMatrix(mat, row, col, op, out)) 824 finally: 825 mat.ops.createsubmatrix = MatCreateSubMatrix_Python 826 return FunctionEnd() 827 cdef Mat sub = None 828 if op == MAT_IGNORE_MATRIX: 829 sub = None 830 elif op == MAT_INITIAL_MATRIX: 831 sub = createSubMatrix(Mat_(mat), IS_(row), IS_(col), None) 832 if sub is not None: 833 addRef(sub.mat) 834 elif op == MAT_REUSE_MATRIX: 835 sub = createSubMatrix(Mat_(mat), IS_(row), IS_(col), Mat_(out[0])) 836 if sub is not None: 837 out[0] = sub.mat 838 return FunctionEnd() 839 840cdef PetscErrorCode MatSetOption_Python( 841 PetscMat mat, 842 PetscMatOption op, 843 PetscBool flag, 844 ) except PETSC_ERR_PYTHON with gil: 845 FunctionBegin(b"MatSetOption_Python") 846 cdef setOption = PyMat(mat).setOption 847 if setOption is not None: 848 setOption(Mat_(mat), <long>op, <bint>(<int>flag)) 849 return FunctionEnd() 850 851cdef PetscErrorCode MatSetUp_Python( 852 PetscMat mat, 853 ) except PETSC_ERR_PYTHON with gil: 854 FunctionBegin(b"MatSetUp_Python") 855 cdef PetscInt rbs = -1, cbs = -1 856 CHKERR(PetscLayoutGetBlockSize(mat.rmap, &rbs)) 857 CHKERR(PetscLayoutGetBlockSize(mat.cmap, &cbs)) 858 if rbs == -1: rbs = 1 859 if cbs == -1: cbs = rbs 860 CHKERR(PetscLayoutSetBlockSize(mat.rmap, rbs)) 861 CHKERR(PetscLayoutSetBlockSize(mat.cmap, cbs)) 862 CHKERR(PetscLayoutSetUp(mat.rmap)) 863 CHKERR(PetscLayoutSetUp(mat.cmap)) 864 mat.preallocated = PETSC_TRUE 865 # 866 cdef char name[2048] 867 cdef PetscBool found = PETSC_FALSE 868 if PyMat(mat).self is None: 869 CHKERR(PetscOptionsGetString(NULL, 870 getPrefix(mat), b"-mat_python_type", 871 name, sizeof(name), &found)) 872 if found and name[0]: 873 CHKERR(MatPythonSetType_PYTHON(mat, name)) 874 if PyMat(mat).self is None: 875 return PetscSETERR(PETSC_ERR_USER, 876 "Python context not set, call one of \n" 877 " * MatPythonSetType(mat, \"[package.]module.class\")\n" 878 " * MatSetFromOptions(mat) and pass option " 879 "-mat_python_type [package.]module.class") 880 # 881 cdef setUp = PyMat(mat).setUp 882 if setUp is not None: 883 setUp(Mat_(mat)) 884 return FunctionEnd() 885 886cdef PetscErrorCode MatAssemblyBegin_Python( 887 PetscMat mat, 888 PetscMatAssemblyType at, 889 ) except PETSC_ERR_PYTHON with gil: 890 FunctionBegin(b"MatAssemblyBegin_Python") 891 cdef assembly = PyMat(mat).assemblyBegin 892 if assembly is not None: 893 assembly(Mat_(mat), <long>at) 894 return FunctionEnd() 895 896cdef PetscErrorCode MatAssemblyEnd_Python( 897 PetscMat mat, 898 PetscMatAssemblyType at, 899 ) except PETSC_ERR_PYTHON with gil: 900 FunctionBegin(b"MatAssemblyEnd_Python") 901 cdef assembly = PyMat(mat).assemblyEnd 902 if assembly is None: 903 assembly = PyMat(mat).assembly 904 if assembly is not None: 905 assembly(Mat_(mat), <long>at) 906 return FunctionEnd() 907 908cdef PetscErrorCode MatZeroEntries_Python( 909 PetscMat mat, 910 ) except PETSC_ERR_PYTHON with gil: 911 FunctionBegin(b"MatZeroEntries_Python") 912 cdef zeroEntries = PyMat(mat).zeroEntries 913 if zeroEntries is None: return UNSUPPORTED(b"zeroEntries") 914 zeroEntries(Mat_(mat)) 915 return FunctionEnd() 916 917cdef PetscErrorCode MatZeroRowsColumns_Python( 918 PetscMat mat, 919 PetscInt numRows, 920 const PetscInt* rows, 921 PetscScalar diag, 922 PetscVec x, 923 PetscVec b, 924 ) except PETSC_ERR_PYTHON with gil: 925 FunctionBegin(b"MatZeroRowsColumns_Python") 926 cdef zeroRowsColumns = PyMat(mat).zeroRowsColumns 927 if zeroRowsColumns is None: return UNSUPPORTED(b"zeroRowsColumns") 928 cdef ndarray pyrows = array_i(numRows, rows) 929 zeroRowsColumns(Mat_(mat), pyrows, toScalar(diag), Vec_(x), Vec_(b)) 930 return FunctionEnd() 931 932cdef PetscErrorCode MatScale_Python( 933 PetscMat mat, 934 PetscScalar s, 935 ) except PETSC_ERR_PYTHON with gil: 936 FunctionBegin(b"MatScale_Python") 937 cdef scale = PyMat(mat).scale 938 if scale is None: return UNSUPPORTED(b"scale") 939 scale(Mat_(mat), toScalar(s)) 940 return FunctionEnd() 941 942cdef PetscErrorCode MatShift_Python( 943 PetscMat mat, 944 PetscScalar s, 945 ) except PETSC_ERR_PYTHON with gil: 946 FunctionBegin(b"MatShift_Python") 947 cdef shift = PyMat(mat).shift 948 if shift is None: return UNSUPPORTED(b"shift") 949 shift(Mat_(mat), toScalar(s)) 950 return FunctionEnd() 951 952cdef PetscErrorCode MatCreateVecs_Python( 953 PetscMat mat, 954 PetscVec *x, 955 PetscVec *y, 956 ) except PETSC_ERR_PYTHON with gil: 957 FunctionBegin(b"MatCreateVecs_Python") 958 cdef createVecs = PyMat(mat).createVecs 959 if createVecs is None: 960 try: 961 mat.ops.getvecs = NULL 962 CHKERR(MatCreateVecs(mat, x, y)) 963 finally: 964 mat.ops.getvecs = MatCreateVecs_Python 965 return FunctionEnd() 966 cdef Vec u, v 967 u, v = createVecs(Mat_(mat)) 968 if x != NULL: 969 x[0] = u.vec 970 u.vec = NULL 971 if y != NULL: 972 y[0] = v.vec 973 v.vec = NULL 974 return FunctionEnd() 975 976cdef PetscErrorCode MatMult_Python( 977 PetscMat mat, 978 PetscVec x, 979 PetscVec y, 980 ) except PETSC_ERR_PYTHON with gil: 981 FunctionBegin(b"MatMult_Python") 982 cdef mult = PyMat(mat).mult 983 if mult is None: return UNSUPPORTED(b"mult") 984 mult(Mat_(mat), Vec_(x), Vec_(y)) 985 return FunctionEnd() 986 987cdef PetscErrorCode MatMultTranspose_Python( 988 PetscMat mat, 989 PetscVec x, 990 PetscVec y, 991 ) except PETSC_ERR_PYTHON with gil: 992 FunctionBegin(b"MatMultTranspose_Python") 993 cdef multTranspose = PyMat(mat).multTranspose 994 if multTranspose is None: 995 try: 996 mat.ops.multtranspose = NULL 997 CHKERR(MatMultTranspose(mat, x, y)) 998 finally: 999 mat.ops.multtranspose = MatMultTranspose_Python 1000 return FunctionEnd() 1001 multTranspose(Mat_(mat), Vec_(x), Vec_(y)) 1002 return FunctionEnd() 1003 1004cdef PetscErrorCode MatMultHermitian_Python( 1005 PetscMat mat, 1006 PetscVec x, 1007 PetscVec y, 1008 ) except PETSC_ERR_PYTHON with gil: 1009 FunctionBegin(b"MatMultHermitian_Python") 1010 cdef multHermitian = PyMat(mat).multHermitian 1011 if multHermitian is None: 1012 try: 1013 mat.ops.multhermitian = NULL 1014 CHKERR(MatMultHermitian(mat, x, y)) 1015 finally: 1016 mat.ops.multhermitian = MatMultHermitian_Python 1017 return FunctionEnd() 1018 multHermitian(Mat_(mat), Vec_(x), Vec_(y)) 1019 return FunctionEnd() 1020 1021cdef PetscErrorCode MatMultAdd_Python( 1022 PetscMat mat, 1023 PetscVec x, 1024 PetscVec v, 1025 PetscVec y, 1026 ) except PETSC_ERR_PYTHON with gil: 1027 FunctionBegin(b"MatMultAdd_Python") 1028 cdef multAdd = PyMat(mat).multAdd 1029 cdef PetscVec t = NULL 1030 if multAdd is None: 1031 if v == y: 1032 CHKERR(VecDuplicate(y, &t)) 1033 CHKERR(MatMult(mat, x, t)) 1034 CHKERR(VecAXPY(y, 1.0, t)) 1035 CHKERR(VecDestroy(&t)) 1036 else: 1037 CHKERR(MatMult(mat, x, y)) 1038 CHKERR(VecAXPY(y, 1.0, v)) 1039 return FunctionEnd() 1040 if multAdd is None: return UNSUPPORTED(b"multAdd") 1041 multAdd(Mat_(mat), Vec_(x), Vec_(v), Vec_(y)) 1042 return FunctionEnd() 1043 1044cdef PetscErrorCode MatMultTransposeAdd_Python( 1045 PetscMat mat, 1046 PetscVec x, 1047 PetscVec v, 1048 PetscVec y, 1049 ) except PETSC_ERR_PYTHON with gil: 1050 FunctionBegin(b"MatMultTransposeAdd_Python") 1051 cdef multTransposeAdd = PyMat(mat).multTransposeAdd 1052 cdef PetscVec t = NULL 1053 if multTransposeAdd is None: 1054 if v == y: 1055 CHKERR(VecDuplicate(y, &t)) 1056 CHKERR(MatMultTranspose(mat, x, t)) 1057 CHKERR(VecAXPY(y, 1.0, t)) 1058 CHKERR(VecDestroy(&t)) 1059 else: 1060 CHKERR(MatMultTranspose(mat, x, y)) 1061 CHKERR(VecAXPY(y, 1.0, v)) 1062 return FunctionEnd() 1063 if multTransposeAdd is None: return UNSUPPORTED(b"multTransposeAdd") 1064 multTransposeAdd(Mat_(mat), Vec_(x), Vec_(v), Vec_(y)) 1065 return FunctionEnd() 1066 1067cdef PetscErrorCode MatMultHermitianAdd_Python( 1068 PetscMat mat, 1069 PetscVec x, 1070 PetscVec v, 1071 PetscVec y, 1072 ) except PETSC_ERR_PYTHON with gil: 1073 FunctionBegin(b"MatMultHermitianAdd_Python") 1074 cdef multHermitianAdd = PyMat(mat).multHermitianAdd 1075 if multHermitianAdd is None: 1076 try: 1077 mat.ops.multhermitianadd = NULL 1078 CHKERR(MatMultHermitianAdd(mat, x, v, y)) 1079 finally: 1080 mat.ops.multhermitianadd = MatMultHermitianAdd_Python 1081 return FunctionEnd() 1082 multHermitianAdd(Mat_(mat), Vec_(x), Vec_(v), Vec_(y)) 1083 return FunctionEnd() 1084 1085cdef PetscErrorCode MatMultDiagonalBlock_Python( 1086 PetscMat mat, 1087 PetscVec x, 1088 PetscVec y, 1089 ) except PETSC_ERR_PYTHON with gil: 1090 FunctionBegin(b"MatMultDiagonalBlock_Python") 1091 cdef multDiagonalBlock = PyMat(mat).multDiagonalBlock 1092 if multDiagonalBlock is None: return UNSUPPORTED(b"multDiagonalBlock") 1093 multDiagonalBlock(Mat_(mat), Vec_(x), Vec_(y)) 1094 return FunctionEnd() 1095 1096cdef PetscErrorCode MatSolve_Python( 1097 PetscMat mat, 1098 PetscVec b, 1099 PetscVec x, 1100 ) except PETSC_ERR_PYTHON with gil: 1101 FunctionBegin(b"MatSolve_Python") 1102 cdef solve = PyMat(mat).solve 1103 if solve is None: return UNSUPPORTED(b"solve") 1104 solve(Mat_(mat), Vec_(b), Vec_(x)) 1105 return FunctionEnd() 1106 1107cdef PetscErrorCode MatSolveTranspose_Python( 1108 PetscMat mat, 1109 PetscVec b, 1110 PetscVec x, 1111 ) except PETSC_ERR_PYTHON with gil: 1112 FunctionBegin(b"MatSolveTranspose_Python") 1113 cdef solveTranspose = PyMat(mat).solveTranspose 1114 if solveTranspose is None: 1115 try: 1116 mat.ops.solvetranspose = NULL 1117 CHKERR(MatSolveTranspose(mat, b, x)) 1118 finally: 1119 mat.ops.solvetranspose = MatSolveTranspose_Python 1120 solveTranspose(Mat_(mat), Vec_(b), Vec_(x)) 1121 return FunctionEnd() 1122 1123cdef PetscErrorCode MatSolveAdd_Python( 1124 PetscMat mat, 1125 PetscVec b, 1126 PetscVec y, 1127 PetscVec x, 1128 ) except PETSC_ERR_PYTHON with gil: 1129 FunctionBegin(b"MatSolveAdd_Python") 1130 cdef solveAdd = PyMat(mat).solveAdd 1131 if solveAdd is None: 1132 try: 1133 mat.ops.solveadd = NULL 1134 CHKERR(MatSolveAdd(mat, b, y, x)) 1135 finally: 1136 mat.ops.solveadd = MatSolveAdd_Python 1137 return FunctionEnd() 1138 solveAdd(Mat_(mat), Vec_(b), Vec_(y), Vec_(x)) 1139 return FunctionEnd() 1140 1141cdef PetscErrorCode MatSolveTransposeAdd_Python( 1142 PetscMat mat, 1143 PetscVec b, 1144 PetscVec y, 1145 PetscVec x, 1146 ) except PETSC_ERR_PYTHON with gil: 1147 FunctionBegin(b"MatSolveTransposeAdd_Python") 1148 cdef solveTransposeAdd = PyMat(mat).solveTransposeAdd 1149 if solveTransposeAdd is None: 1150 try: 1151 mat.ops.solvetransposeadd = NULL 1152 CHKERR(MatSolveTransposeAdd(mat, b, y, x)) 1153 finally: 1154 mat.ops.solvetransposeadd = MatSolveTransposeAdd_Python 1155 return FunctionEnd() 1156 solveTransposeAdd(Mat_(mat), Vec_(b), Vec_(y), Vec_(x)) 1157 return FunctionEnd() 1158 1159cdef PetscErrorCode MatSOR_Python( 1160 PetscMat mat, 1161 PetscVec b, 1162 PetscReal omega, 1163 PetscMatSORType sortype, 1164 PetscReal shift, 1165 PetscInt its, 1166 PetscInt lits, 1167 PetscVec x 1168 ) except PETSC_ERR_PYTHON with gil: 1169 FunctionBegin(b"MatSOR_Python") 1170 cdef SOR = PyMat(mat).SOR 1171 if SOR is None: return UNSUPPORTED(b"SOR") 1172 SOR(Mat_(mat), Vec_(b), asReal(omega), asInt(sortype), asReal(shift), asInt(its), asInt(lits), Vec_(x)) 1173 return FunctionEnd() 1174 1175cdef PetscErrorCode MatGetDiagonal_Python( 1176 PetscMat mat, 1177 PetscVec v, 1178 ) except PETSC_ERR_PYTHON with gil: 1179 FunctionBegin(b"MatGetDiagonal_Python") 1180 cdef getDiagonal = PyMat(mat).getDiagonal 1181 if getDiagonal is None: return UNSUPPORTED(b"getDiagonal") 1182 getDiagonal(Mat_(mat), Vec_(v)) 1183 return FunctionEnd() 1184 1185cdef PetscErrorCode MatSetDiagonal_Python( 1186 PetscMat mat, 1187 PetscVec v, 1188 PetscInsertMode im, 1189 ) except PETSC_ERR_PYTHON with gil: 1190 FunctionBegin(b"MatSetDiagonal_Python") 1191 cdef setDiagonal = PyMat(mat).setDiagonal 1192 cdef bint addv = True if im == PETSC_ADD_VALUES else False 1193 if setDiagonal is None: return UNSUPPORTED(b"setDiagonal") 1194 setDiagonal(Mat_(mat), Vec_(v), addv) 1195 return FunctionEnd() 1196 1197cdef PetscErrorCode MatDiagonalScale_Python( 1198 PetscMat mat, 1199 PetscVec l, 1200 PetscVec r, 1201 ) except PETSC_ERR_PYTHON with gil: 1202 FunctionBegin(b"MatDiagonalScale_Python") 1203 cdef diagonalScale = PyMat(mat).diagonalScale 1204 if diagonalScale is None: return UNSUPPORTED(b"diagonalScale") 1205 diagonalScale(Mat_(mat), Vec_(l), Vec_(r)) 1206 return FunctionEnd() 1207 1208cdef PetscErrorCode MatNorm_Python( 1209 PetscMat mat, 1210 PetscNormType ntype, 1211 PetscReal *nrm, 1212 ) except PETSC_ERR_PYTHON with gil: 1213 FunctionBegin(b"MatNorm_Python") 1214 cdef norm = PyMat(mat).norm 1215 if norm is None: return UNSUPPORTED(b"norm") 1216 retval = norm(Mat_(mat), ntype) 1217 nrm[0] = <PetscReal>retval 1218 return FunctionEnd() 1219 1220cdef PetscErrorCode MatRealPart_Python( 1221 PetscMat mat, 1222 ) except PETSC_ERR_PYTHON with gil: 1223 FunctionBegin(b"MatRealPart_Python") 1224 cdef realPart = PyMat(mat).realPart 1225 if realPart is None: return UNSUPPORTED(b"realPart") 1226 realPart(Mat_(mat)) 1227 return FunctionEnd() 1228 1229cdef PetscErrorCode MatImagPart_Python( 1230 PetscMat mat, 1231 ) except PETSC_ERR_PYTHON with gil: 1232 FunctionBegin(b"MatImagPart_Python") 1233 cdef imagPart = PyMat(mat).imagPart 1234 if imagPart is None: return UNSUPPORTED(b"imagPart") 1235 imagPart(Mat_(mat)) 1236 return FunctionEnd() 1237 1238cdef PetscErrorCode MatConjugate_Python( 1239 PetscMat mat, 1240 ) except PETSC_ERR_PYTHON with gil: 1241 FunctionBegin(b"MatConjugate_Python") 1242 cdef conjugate = PyMat(mat).conjugate 1243 if conjugate is None: return UNSUPPORTED(b"conjugate") 1244 conjugate(Mat_(mat)) 1245 return FunctionEnd() 1246 1247cdef PetscErrorCode MatHasOperation_Python( 1248 PetscMat mat, 1249 PetscMatOperation op, 1250 PetscBool *flag 1251 ) except PETSC_ERR_PYTHON with gil: 1252 FunctionBegin(b"MatHasOperation_Python") 1253 flag[0] = PETSC_FALSE 1254 cdef long i = <long> op 1255 global dMatOps 1256 name = dMatOps.get(i, None) 1257 cdef void** ops = NULL 1258 if name is None: 1259 ops = <void**> mat.ops 1260 if ops and ops[i]: flag[0] = PETSC_TRUE 1261 else: 1262 flag[0] = PETSC_TRUE if getattr(PyMat(mat), name) is not None else PETSC_FALSE 1263 return FunctionEnd() 1264 1265cdef PetscErrorCode MatProductNumeric_Python( 1266 PetscMat mat 1267 ) except PETSC_ERR_PYTHON with gil: 1268 FunctionBegin(b"MatProductNumeric_Python") 1269 cdef PetscMat A = NULL 1270 cdef PetscMat B = NULL 1271 cdef PetscMat C = NULL 1272 cdef PetscMatProductType mtype = MATPRODUCT_UNSPECIFIED 1273 CHKERR(MatProductGetMats(mat, &A, &B, &C)) 1274 CHKERR(MatProductGetType(mat, &mtype)) 1275 1276 mtypes = {MATPRODUCT_AB : 'AB', MATPRODUCT_ABt : 'ABt', MATPRODUCT_AtB : 'AtB', MATPRODUCT_PtAP : 'PtAP', MATPRODUCT_RARt: 'RARt', MATPRODUCT_ABC: 'ABC'} 1277 1278 cdef Mat_Product *product = mat.product 1279 cdef PetscInt i = <PetscInt> <Py_uintptr_t> product.data 1280 if i < 0 or i > 2: 1281 return PetscSETERR(PETSC_ERR_PLIB, 1282 "Corrupted composed id") 1283 cdef PetscMat pM = C if i == 2 else B if i == 1 else A 1284 1285 cdef Mat PyA = Mat_(A) 1286 cdef Mat PyB = Mat_(B) 1287 cdef Mat PyC = Mat_(C) 1288 if mtype == MATPRODUCT_ABC: 1289 mats = (PyA, PyB, PyC) 1290 else: 1291 mats = (PyA, PyB, None) 1292 1293 cdef productNumeric = PyMat(pM).productNumeric 1294 if productNumeric is None: return UNSUPPORTED(b"productNumeric") 1295 productNumeric(PyC if C == pM else PyB if B == pM else PyA, Mat_(mat), mtypes[mtype], *mats) 1296 1297 return FunctionEnd() 1298 1299cdef PetscInt matmatid = -1 1300 1301cdef PetscErrorCode MatProductSymbolic_Python( 1302 PetscMat mat 1303 ) except PETSC_ERR_PYTHON with gil: 1304 FunctionBegin(b"MatProductSymbolic_Python") 1305 cdef PetscMat A = NULL 1306 cdef PetscMat B = NULL 1307 cdef PetscMat C = NULL 1308 cdef PetscMatProductType mtype = MATPRODUCT_UNSPECIFIED 1309 CHKERR(MatProductGetMats(mat, &A, &B, &C)) 1310 CHKERR(MatProductGetType(mat, &mtype)) 1311 1312 mtypes = {MATPRODUCT_AB : 'AB', MATPRODUCT_ABt : 'ABt', MATPRODUCT_AtB : 'AtB', MATPRODUCT_PtAP : 'PtAP', MATPRODUCT_RARt: 'RARt', MATPRODUCT_ABC: 'ABC'} 1313 1314 global matmatid 1315 cdef PetscInt i = -1 1316 cdef PetscBool flg = PETSC_FALSE 1317 CHKERR(PetscObjectComposedDataGetIntPy(<PetscObject>mat, matmatid, &i, &flg)) 1318 if flg is not PETSC_TRUE: 1319 return PetscSETERR(PETSC_ERR_PLIB, 1320 "Missing composed id") 1321 if i < 0 or i > 2: 1322 return PetscSETERR(PETSC_ERR_PLIB, 1323 "Corrupted composed id") 1324 cdef PetscMat pM = C if i == 2 else B if i == 1 else A 1325 1326 cdef Mat PyA = Mat_(A) 1327 cdef Mat PyB = Mat_(B) 1328 cdef Mat PyC = Mat_(C) 1329 if mtype == MATPRODUCT_ABC: 1330 mats = (PyA, PyB, PyC) 1331 else: 1332 mats = (PyA, PyB, None) 1333 1334 cdef productSymbolic = PyMat(pM).productSymbolic 1335 if productSymbolic is None: return UNSUPPORTED(b"productSymbolic") 1336 productSymbolic(PyC if C == pM else PyB if B == pM else PyA, Mat_(mat), mtypes[mtype], *mats) 1337 1338 # Store id in matrix product 1339 cdef Mat_Product *product = mat.product 1340 product.data = <void*> <Py_uintptr_t> i 1341 1342 cdef MatOps ops = mat.ops 1343 ops.productnumeric = MatProductNumeric_Python 1344 return FunctionEnd() 1345 1346cdef PetscErrorCode MatProductSetFromOptions_Python( 1347 PetscMat mat 1348 ) except PETSC_ERR_PYTHON with gil: 1349 FunctionBegin(b"MatProductSetFromOptions_Python") 1350 cdef PetscMat A = NULL 1351 cdef PetscMat B = NULL 1352 cdef PetscMat C = NULL 1353 CHKERR(MatProductGetMats(mat, &A, &B, &C)) 1354 if A == NULL or B == NULL: 1355 return PetscSETERR(PETSC_ERR_PLIB, 1356 "Missing matrices") 1357 1358 cdef PetscMatProductType mtype = MATPRODUCT_UNSPECIFIED 1359 CHKERR(MatProductGetType(mat, &mtype)) 1360 if mtype == MATPRODUCT_UNSPECIFIED: 1361 return PetscSETERR(PETSC_ERR_PLIB, 1362 "Unknown product type") 1363 1364 mtypes = {MATPRODUCT_AB : 'AB', MATPRODUCT_ABt : 'ABt', MATPRODUCT_AtB : 'AtB', MATPRODUCT_PtAP : 'PtAP', MATPRODUCT_RARt: 'RARt', MATPRODUCT_ABC: 'ABC'} 1365 1366 cdef Mat PyA = Mat_(A) 1367 cdef Mat PyB = Mat_(B) 1368 cdef Mat PyC = Mat_(C) 1369 if mtype == MATPRODUCT_ABC: 1370 mats = (PyA, PyB, PyC) 1371 else: 1372 mats = (PyA, PyB, None) 1373 1374 # Find Python matrix in mats able to perform the product 1375 found = False 1376 cdef PetscBool mispy = PETSC_FALSE 1377 cdef PetscMat pM = NULL 1378 cdef Mat mm 1379 cdef PetscInt i = -1 1380 for i in range(len(mats)): 1381 if mats[i] is None: continue 1382 mm = mats[i] 1383 pM = <PetscMat>mm.mat 1384 CHKERR(PetscObjectTypeCompare(<PetscObject>pM, MATPYTHON, &mispy)) 1385 if mispy: 1386 if PyMat(pM).productSetFromOptions is not None: 1387 found = PyMat(pM).productSetFromOptions(PyC if C == pM else PyB if B == pM else PyA, mtypes[mtype], *mats) 1388 if found: break 1389 if not found: 1390 return FunctionEnd() 1391 1392 cdef MatOps ops = mat.ops 1393 ops.productsymbolic = MatProductSymbolic_Python 1394 1395 # Store index (within the product) of the Python matrix which is capable of performing the operation 1396 # Cannot be stored in mat.product.data at this stage 1397 # Symbolic operation will get this index and store it in the product data 1398 global matmatid 1399 if matmatid < 0: 1400 CHKERR(PetscObjectComposedDataRegisterPy(&matmatid)) 1401 CHKERR(PetscObjectComposedDataSetIntPy(<PetscObject>mat, matmatid, i)) 1402 1403 return FunctionEnd() 1404 1405# -------------------------------------------------------------------- 1406 1407cdef extern from * nogil: 1408 struct _PCOps: 1409 PetscErrorCode (*destroy)(PetscPC) except PETSC_ERR_PYTHON 1410 PetscErrorCode (*setup)(PetscPC) except PETSC_ERR_PYTHON 1411 PetscErrorCode (*reset)(PetscPC) except PETSC_ERR_PYTHON 1412 PetscErrorCode (*setfromoptions)(PetscPC, PetscOptionItems) except PETSC_ERR_PYTHON 1413 PetscErrorCode (*view)(PetscPC, PetscViewer) except PETSC_ERR_PYTHON 1414 PetscErrorCode (*presolve)(PetscPC, PetscKSP, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1415 PetscErrorCode (*postsolve)(PetscPC, PetscKSP, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1416 PetscErrorCode (*apply)(PetscPC, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1417 PetscErrorCode (*matapply)(PetscPC, PetscMat, PetscMat) except PETSC_ERR_PYTHON 1418 PetscErrorCode (*applytranspose)(PetscPC, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1419 PetscErrorCode (*matapplytranspose)(PetscPC, PetscMat, PetscMat) except PETSC_ERR_PYTHON 1420 PetscErrorCode (*applysymmetricleft)(PetscPC, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1421 PetscErrorCode (*applysymmetricright)(PetscPC, PetscVec, PetscVec) except PETSC_ERR_PYTHON 1422 ctypedef _PCOps *PCOps 1423 struct _p_PC: 1424 void *data 1425 PCOps ops 1426 1427 1428@cython.internal 1429cdef class _PyPC(_PyObj): pass 1430cdef inline _PyPC PyPC(PetscPC pc): 1431 if pc != NULL and pc.data != NULL: 1432 return <_PyPC>pc.data 1433 else: 1434 return _PyPC.__new__(_PyPC) 1435 1436cdef public PetscErrorCode PCPythonGetContext(PetscPC pc, void **ctx) \ 1437 except PETSC_ERR_PYTHON: 1438 FunctionBegin(b"PCPythonGetContext") 1439 PyPC(pc).getcontext(ctx) 1440 return FunctionEnd() 1441 1442cdef public PetscErrorCode PCPythonSetContext(PetscPC pc, void *ctx) \ 1443 except PETSC_ERR_PYTHON: 1444 FunctionBegin(b"PCPythonSetContext") 1445 PyPC(pc).setcontext(ctx, PC_(pc)) 1446 return FunctionEnd() 1447 1448cdef PetscErrorCode PCPythonSetType_PYTHON(PetscPC pc, const char *name) \ 1449 except PETSC_ERR_PYTHON with gil: 1450 FunctionBegin(b"PCPythonSetType_PYTHON") 1451 if name == NULL: return FunctionEnd() # XXX 1452 cdef object ctx = createcontext(name) 1453 PCPythonSetContext(pc, <void*>ctx) 1454 PyPC(pc).setname(name) 1455 return FunctionEnd() 1456 1457cdef PetscErrorCode PCPythonGetType_PYTHON(PetscPC pc, const char *name[]) \ 1458 except PETSC_ERR_PYTHON with gil: 1459 FunctionBegin(b"PCPythonGetType_PYTHON") 1460 name[0] = PyPC(pc).getname() 1461 return FunctionEnd() 1462 1463cdef PetscErrorCode PCCreate_Python( 1464 PetscPC pc, 1465 ) except PETSC_ERR_PYTHON with gil: 1466 FunctionBegin(b"PCCreate_Python") 1467 cdef PCOps ops = pc.ops 1468 ops.reset = PCReset_Python 1469 ops.destroy = PCDestroy_Python 1470 ops.setup = PCSetUp_Python 1471 ops.setfromoptions = PCSetFromOptions_Python 1472 ops.view = PCView_Python 1473 ops.presolve = PCPreSolve_Python 1474 ops.postsolve = PCPostSolve_Python 1475 ops.apply = PCApply_Python 1476 ops.matapply = PCMatApply_Python 1477 ops.applytranspose = PCApplyTranspose_Python 1478 ops.matapplytranspose = PCMatApplyTranspose_Python 1479 ops.applysymmetricleft = PCApplySymmetricLeft_Python 1480 ops.applysymmetricright = PCApplySymmetricRight_Python 1481 # 1482 CHKERR(PetscObjectComposeFunction( 1483 <PetscObject>pc, b"PCPythonSetType_C", 1484 <PetscVoidFunction>PCPythonSetType_PYTHON)) 1485 CHKERR(PetscObjectComposeFunction( 1486 <PetscObject>pc, b"PCPythonGetType_C", 1487 <PetscVoidFunction>PCPythonGetType_PYTHON)) 1488 # 1489 cdef ctx = PyPC(NULL) 1490 pc.data = <void*> ctx 1491 Py_INCREF(<PyObject*>pc.data) 1492 return FunctionEnd() 1493 1494cdef inline PetscErrorCode PCDestroy_Python_inner( 1495 PetscPC pc, 1496 ) except PETSC_ERR_PYTHON with gil: 1497 try: 1498 addRef(pc) 1499 PCPythonSetContext(pc, NULL) 1500 finally: 1501 delRef(pc) 1502 Py_DECREF(<PyObject*>pc.data) 1503 pc.data = NULL 1504 return PETSC_SUCCESS 1505 1506cdef PetscErrorCode PCDestroy_Python( 1507 PetscPC pc, 1508 ) except PETSC_ERR_PYTHON nogil: 1509 FunctionBegin(b"PCDestroy_Python") 1510 CHKERR(PetscObjectComposeFunction( 1511 <PetscObject>pc, b"PCPythonSetType_C", 1512 <PetscVoidFunction>NULL)) 1513 CHKERR(PetscObjectComposeFunction( 1514 <PetscObject>pc, b"PCPythonGetType_C", 1515 <PetscVoidFunction>NULL)) 1516 # 1517 if Py_IsInitialized(): PCDestroy_Python_inner(pc) 1518 return FunctionEnd() 1519 1520cdef PetscErrorCode PCSetUp_Python( 1521 PetscPC pc, 1522 ) except PETSC_ERR_PYTHON with gil: 1523 FunctionBegin(b"PCSetUp_Python") 1524 cdef char name[2048] 1525 cdef PetscBool found = PETSC_FALSE 1526 if PyPC(pc).self is None: 1527 CHKERR(PetscOptionsGetString(NULL, 1528 getPrefix(pc), b"-pc_python_type", 1529 name, sizeof(name), &found)) 1530 if found and name[0]: 1531 CHKERR(PCPythonSetType_PYTHON(pc, name)) 1532 if PyPC(pc).self is None: 1533 return PetscSETERR(PETSC_ERR_USER, 1534 "Python context not set, call one of \n" 1535 " * PCPythonSetType(pc, \"[package.]module.class\")\n" 1536 " * PCSetFromOptions(pc) and pass option " 1537 "-pc_python_type [package.]module.class") 1538 # 1539 cdef setUp = PyPC(pc).setUp 1540 if setUp is not None: 1541 setUp(PC_(pc)) 1542 # 1543 cdef o = PyPC(pc) 1544 cdef PCOps ops = pc.ops 1545 if o.applyTranspose is None: 1546 ops.applytranspose = NULL 1547 if o.applySymmetricLeft is None: 1548 ops.applysymmetricleft = NULL 1549 if o.applySymmetricRight is None: 1550 ops.applysymmetricright = NULL 1551 # 1552 return FunctionEnd() 1553 1554cdef inline PetscErrorCode PCReset_Python_inner( 1555 PetscPC pc, 1556 ) except PETSC_ERR_PYTHON with gil: 1557 cdef reset = PyPC(pc).reset 1558 if reset is not None: 1559 reset(PC_(pc)) 1560 return PETSC_SUCCESS 1561 1562cdef PetscErrorCode PCReset_Python( 1563 PetscPC pc, 1564 ) except PETSC_ERR_PYTHON nogil: 1565 if getRef(pc) == 0: return PETSC_SUCCESS 1566 FunctionBegin(b"PCReset_Python") 1567 if Py_IsInitialized(): PCReset_Python_inner(pc) 1568 return FunctionEnd() 1569 1570cdef PetscErrorCode PCSetFromOptions_Python( 1571 PetscPC pc, 1572 PetscOptionItems PetscOptionsObject, 1573 ) except PETSC_ERR_PYTHON with gil: 1574 FunctionBegin(b"PCSetFromOptions_Python") 1575 cdef char name[2048], *defval = PyPC(pc).getname() 1576 cdef PetscBool found = PETSC_FALSE 1577 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 1578 CHKERR(PetscOptionsString( 1579 b"-pc_python_type", b"Python [package.]module[.{class|function}]", 1580 b"PCPythonSetType", defval, name, sizeof(name), &found)); <void>opts 1581 if found and name[0]: 1582 CHKERR(PCPythonSetType_PYTHON(pc, name)) 1583 # 1584 cdef setFromOptions = PyPC(pc).setFromOptions 1585 if setFromOptions is not None: 1586 setFromOptions(PC_(pc)) 1587 return FunctionEnd() 1588 1589cdef PetscErrorCode PCView_Python( 1590 PetscPC pc, 1591 PetscViewer vwr, 1592 ) except PETSC_ERR_PYTHON with gil: 1593 FunctionBegin(b"PCView_Python") 1594 viewcontext(PyPC(pc), vwr) 1595 cdef view = PyPC(pc).view 1596 if view is not None: 1597 view(PC_(pc), Viewer_(vwr)) 1598 return FunctionEnd() 1599 1600cdef PetscErrorCode PCPreSolve_Python( 1601 PetscPC pc, 1602 PetscKSP ksp, 1603 PetscVec b, 1604 PetscVec x, 1605 ) except PETSC_ERR_PYTHON with gil: 1606 FunctionBegin(b"PCPreSolve_Python") 1607 cdef preSolve = PyPC(pc).preSolve 1608 if preSolve is not None: 1609 preSolve(PC_(pc), KSP_(ksp), Vec_(b), Vec_(x)) 1610 return FunctionEnd() 1611 1612cdef PetscErrorCode PCPostSolve_Python( 1613 PetscPC pc, 1614 PetscKSP ksp, 1615 PetscVec b, 1616 PetscVec x, 1617 ) except PETSC_ERR_PYTHON with gil: 1618 FunctionBegin(b"PCPostSolve_Python") 1619 cdef postSolve = PyPC(pc).postSolve 1620 if postSolve is not None: 1621 postSolve(PC_(pc), KSP_(ksp), Vec_(b), Vec_(x)) 1622 return FunctionEnd() 1623 1624cdef PetscErrorCode PCApply_Python( 1625 PetscPC pc, 1626 PetscVec x, 1627 PetscVec y, 1628 ) except PETSC_ERR_PYTHON with gil: 1629 FunctionBegin(b"PCApply_Python") 1630 cdef apply = PyPC(pc).apply 1631 apply(PC_(pc), Vec_(x), Vec_(y)) 1632 return FunctionEnd() 1633 1634cdef PetscErrorCode PCApplyTranspose_Python( 1635 PetscPC pc, 1636 PetscVec x, 1637 PetscVec y, 1638 ) except PETSC_ERR_PYTHON with gil: 1639 FunctionBegin(b"PCApplyTranspose_Python") 1640 cdef applyTranspose = PyPC(pc).applyTranspose 1641 applyTranspose(PC_(pc), Vec_(x), Vec_(y)) 1642 return FunctionEnd() 1643 1644cdef PetscErrorCode PCApplySymmetricLeft_Python( 1645 PetscPC pc, 1646 PetscVec x, 1647 PetscVec y, 1648 ) except PETSC_ERR_PYTHON with gil: 1649 FunctionBegin(b"PCApplySymmetricLeft_Python") 1650 cdef applySymmetricLeft = PyPC(pc).applySymmetricLeft 1651 applySymmetricLeft(PC_(pc), Vec_(x), Vec_(y)) 1652 return FunctionEnd() 1653 1654cdef PetscErrorCode PCApplySymmetricRight_Python( 1655 PetscPC pc, 1656 PetscVec x, 1657 PetscVec y, 1658 ) except PETSC_ERR_PYTHON with gil: 1659 FunctionBegin(b"PCApplySymmetricRight_Python") 1660 cdef applySymmetricRight = PyPC(pc).applySymmetricRight 1661 applySymmetricRight(PC_(pc), Vec_(x), Vec_(y)) 1662 return FunctionEnd() 1663 1664cdef PetscErrorCode PCMatApply_Python( 1665 PetscPC pc, 1666 PetscMat X, 1667 PetscMat Y, 1668 ) except PETSC_ERR_PYTHON with gil: 1669 FunctionBegin(b"PCMatApply_Python") 1670 cdef matApply = PyPC(pc).matApply 1671 if matApply is None: 1672 try: 1673 pc.ops.matapply = NULL 1674 CHKERR(PCMatApply(pc, X, Y)) 1675 finally: 1676 pc.ops.matapply = PCMatApply_Python 1677 return FunctionEnd() 1678 matApply(PC_(pc), Mat_(X), Mat_(Y)) 1679 return FunctionEnd() 1680 1681cdef PetscErrorCode PCMatApplyTranspose_Python( 1682 PetscPC pc, 1683 PetscMat X, 1684 PetscMat Y, 1685 ) except PETSC_ERR_PYTHON with gil: 1686 FunctionBegin(b"PCMatApplyTranspose_Python") 1687 cdef matApplyTranspose = PyPC(pc).matApplyTranspose 1688 if matApplyTranspose is None: 1689 try: 1690 pc.ops.matapplytranspose = NULL 1691 CHKERR(PCMatApplyTranspose(pc, X, Y)) 1692 finally: 1693 pc.ops.matapplytranspose = PCMatApplyTranspose_Python 1694 return FunctionEnd() 1695 matApplyTranspose(PC_(pc), Mat_(X), Mat_(Y)) 1696 return FunctionEnd() 1697 1698# -------------------------------------------------------------------- 1699 1700cdef extern from * nogil: 1701 struct _KSPOps: 1702 PetscErrorCode (*destroy)(PetscKSP) except PETSC_ERR_PYTHON 1703 PetscErrorCode (*setup)(PetscKSP) except PETSC_ERR_PYTHON 1704 PetscErrorCode (*reset)(PetscKSP) except PETSC_ERR_PYTHON 1705 PetscErrorCode (*setfromoptions)(PetscKSP, PetscOptionItems) except PETSC_ERR_PYTHON 1706 PetscErrorCode (*view)(PetscKSP, PetscViewer) except PETSC_ERR_PYTHON 1707 PetscErrorCode (*solve)(PetscKSP) except PETSC_ERR_PYTHON 1708 PetscErrorCode (*buildsolution)(PetscKSP, PetscVec, PetscVec*) except PETSC_ERR_PYTHON 1709 PetscErrorCode (*buildresidual)(PetscKSP, PetscVec, PetscVec, PetscVec*) except PETSC_ERR_PYTHON 1710 ctypedef _KSPOps *KSPOps 1711 struct _p_KSP: 1712 void *data 1713 KSPOps ops 1714 PetscBool transpose_solve 1715 PetscInt iter"its", max_its"max_it" 1716 PetscReal norm"rnorm" 1717 PetscKSPConvergedReason reason 1718 1719cdef extern from * nogil: # custom.h 1720 PetscErrorCode KSPConverged(PetscKSP, PetscInt, PetscReal, PetscKSPConvergedReason*) 1721 PetscErrorCode KSPLogHistory(PetscKSP, PetscReal) 1722 1723 1724@cython.internal 1725cdef class _PyKSP(_PyObj): pass 1726cdef inline _PyKSP PyKSP(PetscKSP ksp): 1727 if ksp != NULL and ksp.data != NULL: 1728 return <_PyKSP>ksp.data 1729 else: 1730 return _PyKSP.__new__(_PyKSP) 1731 1732cdef public PetscErrorCode KSPPythonGetContext(PetscKSP ksp, void **ctx) \ 1733 except PETSC_ERR_PYTHON: 1734 FunctionBegin(b"KSPPythonGetContext") 1735 PyKSP(ksp).getcontext(ctx) 1736 return FunctionEnd() 1737 1738cdef public PetscErrorCode KSPPythonSetContext(PetscKSP ksp, void *ctx) \ 1739 except PETSC_ERR_PYTHON: 1740 FunctionBegin(b"KSPPythonSetContext") 1741 PyKSP(ksp).setcontext(ctx, KSP_(ksp)) 1742 return FunctionEnd() 1743 1744cdef PetscErrorCode KSPPythonSetType_PYTHON(PetscKSP ksp, const char *name) \ 1745 except PETSC_ERR_PYTHON with gil: 1746 FunctionBegin(b"KSPPythonSetType_PYTHON") 1747 if name == NULL: return FunctionEnd() # XXX 1748 cdef object ctx = createcontext(name) 1749 KSPPythonSetContext(ksp, <void*>ctx) 1750 PyKSP(ksp).setname(name) 1751 return FunctionEnd() 1752 1753cdef PetscErrorCode KSPPythonGetType_PYTHON(PetscKSP ksp, const char *name[]) \ 1754 except PETSC_ERR_PYTHON with gil: 1755 FunctionBegin(b"KSPPythonGetType_PYTHON") 1756 name[0] = PyKSP(ksp).getname() 1757 return FunctionEnd() 1758 1759cdef PetscErrorCode KSPCreate_Python( 1760 PetscKSP ksp, 1761 ) except PETSC_ERR_PYTHON with gil: 1762 FunctionBegin(b"KSPCreate_Python") 1763 cdef KSPOps ops = ksp.ops 1764 ops.reset = KSPReset_Python 1765 ops.destroy = KSPDestroy_Python 1766 ops.setup = KSPSetUp_Python 1767 ops.setfromoptions = KSPSetFromOptions_Python 1768 ops.view = KSPView_Python 1769 ops.solve = KSPSolve_Python 1770 ops.buildsolution = KSPBuildSolution_Python 1771 ops.buildresidual = KSPBuildResidual_Python 1772 # 1773 CHKERR(PetscObjectComposeFunction( 1774 <PetscObject>ksp, b"KSPPythonSetType_C", 1775 <PetscVoidFunction>KSPPythonSetType_PYTHON)) 1776 CHKERR(PetscObjectComposeFunction( 1777 <PetscObject>ksp, b"KSPPythonGetType_C", 1778 <PetscVoidFunction>KSPPythonGetType_PYTHON)) 1779 # 1780 cdef ctx = PyKSP(NULL) 1781 ksp.data = <void*> ctx 1782 Py_INCREF(<PyObject*>ksp.data) 1783 # 1784 CHKERR(KSPSetSupportedNorm( 1785 ksp, KSP_NORM_PRECONDITIONED, PC_LEFT, 3)) 1786 CHKERR(KSPSetSupportedNorm( 1787 ksp, KSP_NORM_UNPRECONDITIONED, PC_RIGHT, 3)) 1788 CHKERR(KSPSetSupportedNorm( 1789 ksp, KSP_NORM_UNPRECONDITIONED, PC_LEFT, 2)) 1790 CHKERR(KSPSetSupportedNorm( 1791 ksp, KSP_NORM_PRECONDITIONED, PC_RIGHT, 2)) 1792 CHKERR(KSPSetSupportedNorm( 1793 ksp, KSP_NORM_PRECONDITIONED, PC_SYMMETRIC, 1)) 1794 CHKERR(KSPSetSupportedNorm( 1795 ksp, KSP_NORM_UNPRECONDITIONED, PC_SYMMETRIC, 1)) 1796 return FunctionEnd() 1797 1798cdef inline PetscErrorCode KSPDestroy_Python_inner( 1799 PetscKSP ksp, 1800 ) except PETSC_ERR_PYTHON with gil: 1801 try: 1802 addRef(ksp) 1803 KSPPythonSetContext(ksp, NULL) 1804 finally: 1805 delRef(ksp) 1806 Py_DECREF(<PyObject*>ksp.data) 1807 ksp.data = NULL 1808 return PETSC_SUCCESS 1809 1810cdef PetscErrorCode KSPDestroy_Python( 1811 PetscKSP ksp, 1812 ) except PETSC_ERR_PYTHON nogil: 1813 FunctionBegin(b"KSPDestroy_Python") 1814 CHKERR(PetscObjectComposeFunction( 1815 <PetscObject>ksp, b"KSPPythonSetType_C", 1816 <PetscVoidFunction>NULL)) 1817 CHKERR(PetscObjectComposeFunction( 1818 <PetscObject>ksp, b"KSPPythonGetType_C", 1819 <PetscVoidFunction>NULL)) 1820 # 1821 if Py_IsInitialized(): KSPDestroy_Python_inner(ksp) 1822 return FunctionEnd() 1823 1824cdef PetscErrorCode KSPSetUp_Python( 1825 PetscKSP ksp, 1826 ) except PETSC_ERR_PYTHON with gil: 1827 FunctionBegin(b"KSPSetUp_Python") 1828 cdef char name[2048] 1829 cdef PetscBool found = PETSC_FALSE 1830 if PyKSP(ksp).self is None: 1831 CHKERR(PetscOptionsGetString(NULL, 1832 getPrefix(ksp), b"-ksp_python_type", 1833 name, sizeof(name), &found)) 1834 if found and name[0]: 1835 CHKERR(KSPPythonSetType_PYTHON(ksp, name)) 1836 if PyKSP(ksp).self is None: 1837 return PetscSETERR(PETSC_ERR_USER, 1838 "Python context not set, call one of \n" 1839 " * KSPPythonSetType(ksp, \"[package.]module.class\")\n" 1840 " * KSPSetFromOptions(ksp) and pass option " 1841 "-ksp_python_type [package.]module.class") 1842 # 1843 cdef setUp = PyKSP(ksp).setUp 1844 if setUp is not None: 1845 setUp(KSP_(ksp)) 1846 return FunctionEnd() 1847 1848cdef inline PetscErrorCode KSPReset_Python_inner( 1849 PetscKSP ksp, 1850 ) except PETSC_ERR_PYTHON with gil: 1851 cdef reset = PyKSP(ksp).reset 1852 if reset is not None: 1853 reset(KSP_(ksp)) 1854 return PETSC_SUCCESS 1855 1856cdef PetscErrorCode KSPReset_Python( 1857 PetscKSP ksp, 1858 ) except PETSC_ERR_PYTHON nogil: 1859 if getRef(ksp) == 0: return PETSC_SUCCESS 1860 FunctionBegin(b"KSPReset_Python") 1861 CHKERR(PetscObjectCompose(<PetscObject>ksp, b"@ksp.vec_work_sol", NULL)) 1862 CHKERR(PetscObjectCompose(<PetscObject>ksp, b"@ksp.vec_work_res", NULL)) 1863 if Py_IsInitialized(): KSPReset_Python_inner(ksp) 1864 return FunctionEnd() 1865 1866cdef PetscErrorCode KSPSetFromOptions_Python( 1867 PetscKSP ksp, 1868 PetscOptionItems PetscOptionsObject 1869 ) except PETSC_ERR_PYTHON with gil: 1870 FunctionBegin(b"KSPSetFromOptions_Python") 1871 cdef char name[2048], *defval = PyKSP(ksp).getname() 1872 cdef PetscBool found = PETSC_FALSE 1873 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 1874 CHKERR(PetscOptionsString( 1875 b"-ksp_python_type", b"Python [package.]module[.{class|function}]", 1876 b"KSPPythonSetType", defval, name, sizeof(name), &found)); <void>opts 1877 if found and name[0]: 1878 CHKERR(KSPPythonSetType_PYTHON(ksp, name)) 1879 # 1880 cdef setFromOptions = PyKSP(ksp).setFromOptions 1881 if setFromOptions is not None: 1882 setFromOptions(KSP_(ksp)) 1883 return FunctionEnd() 1884 1885cdef PetscErrorCode KSPView_Python( 1886 PetscKSP ksp, 1887 PetscViewer vwr, 1888 ) except PETSC_ERR_PYTHON with gil: 1889 FunctionBegin(b"KSPView_Python") 1890 viewcontext(PyKSP(ksp), vwr) 1891 cdef view = PyKSP(ksp).view 1892 if view is not None: 1893 view(KSP_(ksp), Viewer_(vwr)) 1894 return FunctionEnd() 1895 1896cdef PetscErrorCode KSPBuildSolution_Python( 1897 PetscKSP ksp, 1898 PetscVec v, 1899 PetscVec *V, 1900 ) except PETSC_ERR_PYTHON with gil: 1901 FunctionBegin(b"KSPBuildSolution_Python") 1902 cdef PetscVec x = v 1903 cdef buildSolution = PyKSP(ksp).buildSolution 1904 if buildSolution is not None: 1905 if x == NULL: pass # XXX 1906 buildSolution(KSP_(ksp), Vec_(x)) 1907 if V != NULL: V[0] = x 1908 else: 1909 CHKERR(KSPBuildSolutionDefault(ksp, v, V)) 1910 return FunctionEnd() 1911 1912cdef PetscErrorCode KSPBuildResidual_Python( 1913 PetscKSP ksp, 1914 PetscVec t, 1915 PetscVec v, 1916 PetscVec *V, 1917 ) except PETSC_ERR_PYTHON with gil: 1918 FunctionBegin(b"KSPBuildResidual_Python") 1919 cdef buildResidual = PyKSP(ksp).buildResidual 1920 if buildResidual is not None: 1921 buildResidual(KSP_(ksp), Vec_(t), Vec_(v)) 1922 if V != NULL: V[0] = v 1923 else: 1924 CHKERR(KSPBuildResidualDefault(ksp, t, v, V)) 1925 return FunctionEnd() 1926 1927cdef PetscErrorCode KSPSolve_Python( 1928 PetscKSP ksp, 1929 ) except PETSC_ERR_PYTHON with gil: 1930 FunctionBegin(b"KSPSolve_Python") 1931 cdef PetscVec B = NULL, X = NULL 1932 CHKERR(KSPGetRhs(ksp, &B)) 1933 CHKERR(KSPGetSolution(ksp, &X)) 1934 # 1935 ksp.iter = 0 1936 ksp.reason = KSP_CONVERGED_ITERATING 1937 # 1938 cdef solve = None 1939 if ksp.transpose_solve: 1940 solve = PyKSP(ksp).solveTranspose 1941 else: 1942 solve = PyKSP(ksp).solve 1943 if solve is not None: 1944 solve(KSP_(ksp), Vec_(B), Vec_(X)) 1945 else: 1946 KSPSolve_Python_default(ksp, B, X) 1947 return FunctionEnd() 1948 1949cdef PetscErrorCode KSPSolve_Python_default( 1950 PetscKSP ksp, 1951 PetscVec B, 1952 PetscVec X, 1953 ) except PETSC_ERR_PYTHON with gil: 1954 FunctionBegin(b"KSPSolve_Python_default") 1955 cdef PetscVec t = NULL 1956 CHKERR(PetscObjectQuery( 1957 <PetscObject>ksp, 1958 b"@ksp.vec_work_sol", 1959 <PetscObject*>&t)) 1960 if t == NULL: 1961 CHKERR(VecDuplicate(X, &t)) 1962 CHKERR(PetscObjectCompose( 1963 <PetscObject>ksp, 1964 b"@ksp.vec_work_sol", 1965 <PetscObject>t)) 1966 cdef PetscVec v = NULL 1967 CHKERR(PetscObjectQuery( 1968 <PetscObject>ksp, 1969 b"@ksp.vec_work_res", 1970 <PetscObject*>&v)) 1971 if v == NULL: 1972 CHKERR(VecDuplicate(B, &v)) 1973 CHKERR(PetscObjectCompose( 1974 <PetscObject>ksp, 1975 b"@ksp.vec_work_res", 1976 <PetscObject>v)) 1977 # 1978 cdef PetscVec R = NULL 1979 cdef PetscReal rnorm = 0 1980 # 1981 CHKERR(KSPBuildResidual(ksp, t, v, &R)) 1982 CHKERR(VecNorm(R, PETSC_NORM_2, &rnorm)) 1983 # 1984 CHKERR(KSPConverged(ksp, ksp.iter, rnorm, &ksp.reason)) 1985 CHKERR(KSPLogHistory(ksp, ksp.norm)) 1986 CHKERR(KSPMonitor(ksp, ksp.iter, ksp.norm)) 1987 for its from 0 <= its < ksp.max_its: 1988 <void> its # unused 1989 if ksp.reason: break 1990 KSPPreStep_Python(ksp) 1991 # 1992 KSPStep_Python(ksp, B, X) # FIXME? B? 1993 CHKERR(KSPBuildResidual(ksp, t, v, &R)) 1994 CHKERR(VecNorm(R, PETSC_NORM_2, &rnorm)) 1995 ksp.iter += 1 1996 # 1997 KSPPostStep_Python(ksp) 1998 CHKERR(KSPConverged(ksp, ksp.iter, rnorm, &ksp.reason)) 1999 CHKERR(KSPLogHistory(ksp, ksp.norm)) 2000 CHKERR(KSPMonitor(ksp, ksp.iter, ksp.norm)) 2001 if ksp.iter == ksp.max_its: 2002 if ksp.reason == KSP_CONVERGED_ITERATING: 2003 ksp.reason = KSP_DIVERGED_MAX_IT 2004 # 2005 return FunctionEnd() 2006 2007cdef PetscErrorCode KSPPreStep_Python( 2008 PetscKSP ksp, 2009 ) except PETSC_ERR_PYTHON with gil: 2010 FunctionBegin(b"KSPPreStep_Python") 2011 cdef preStep = PyKSP(ksp).preStep 2012 if preStep is not None: 2013 preStep(KSP_(ksp)) 2014 return FunctionEnd() 2015 2016cdef PetscErrorCode KSPPostStep_Python( 2017 PetscKSP ksp, 2018 ) except PETSC_ERR_PYTHON with gil: 2019 FunctionBegin(b"KSPPostStep_Python") 2020 cdef postStep = PyKSP(ksp).postStep 2021 if postStep is not None: 2022 postStep(KSP_(ksp)) 2023 return FunctionEnd() 2024 2025cdef PetscErrorCode KSPStep_Python( 2026 PetscKSP ksp, 2027 PetscVec B, 2028 PetscVec X, 2029 ) except PETSC_ERR_PYTHON with gil: 2030 FunctionBegin(b"KSPStep_Python") 2031 cdef step = None 2032 if ksp.transpose_solve: 2033 step = PyKSP(ksp).stepTranspose 2034 if step is None: return UNSUPPORTED(b"stepTranspose") 2035 else: 2036 step = PyKSP(ksp).step 2037 if step is None: return UNSUPPORTED(b"step") 2038 step(KSP_(ksp), Vec_(B), Vec_(X)) 2039 return FunctionEnd() 2040 2041# -------------------------------------------------------------------- 2042 2043cdef extern from * nogil: 2044 struct _SNESOps: 2045 PetscErrorCode (*destroy)(PetscSNES) except PETSC_ERR_PYTHON 2046 PetscErrorCode (*setup)(PetscSNES) except PETSC_ERR_PYTHON 2047 PetscErrorCode (*reset)(PetscSNES) except PETSC_ERR_PYTHON 2048 PetscErrorCode (*setfromoptions)(PetscSNES, PetscOptionItems) except PETSC_ERR_PYTHON 2049 PetscErrorCode (*view)(PetscSNES, PetscViewer) except PETSC_ERR_PYTHON 2050 PetscErrorCode (*solve)(PetscSNES) except PETSC_ERR_PYTHON 2051 ctypedef _SNESOps *SNESOps 2052 struct _p_SNES: 2053 void *data 2054 SNESOps ops 2055 PetscInt iter, max_its, linear_its 2056 PetscReal norm, xnorm, ynorm, rtol, ttol 2057 PetscSNESConvergedReason reason 2058 PetscVec vec_sol, vec_sol_update, vec_func 2059 PetscMat jacobian, jacobian_pre 2060 PetscKSP ksp 2061 2062cdef extern from * nogil: # custom.h 2063 PetscErrorCode SNESLogHistory(PetscSNES, PetscReal, PetscInt) 2064 PetscErrorCode SNESComputeUpdate(PetscSNES) 2065 2066 2067@cython.internal 2068cdef class _PySNES(_PyObj): pass 2069cdef inline _PySNES PySNES(PetscSNES snes): 2070 if snes != NULL and snes.data != NULL: 2071 return <_PySNES>snes.data 2072 else: 2073 return _PySNES.__new__(_PySNES) 2074 2075cdef public PetscErrorCode SNESPythonGetContext(PetscSNES snes, void **ctx) \ 2076 except PETSC_ERR_PYTHON: 2077 FunctionBegin(b"SNESPythonGetContext ") 2078 PySNES(snes).getcontext(ctx) 2079 return FunctionEnd() 2080 2081cdef public PetscErrorCode SNESPythonSetContext(PetscSNES snes, void *ctx) \ 2082 except PETSC_ERR_PYTHON: 2083 FunctionBegin(b"SNESPythonSetContext ") 2084 PySNES(snes).setcontext(ctx, SNES_(snes)) 2085 return FunctionEnd() 2086 2087cdef PetscErrorCode SNESPythonSetType_PYTHON(PetscSNES snes, const char *name) \ 2088 except PETSC_ERR_PYTHON with gil: 2089 FunctionBegin(b"SNESPythonSetType_PYTHON") 2090 if name == NULL: return FunctionEnd() # XXX 2091 cdef object ctx = createcontext(name) 2092 SNESPythonSetContext(snes, <void*>ctx) 2093 PySNES(snes).setname(name) 2094 return FunctionEnd() 2095 2096cdef PetscErrorCode SNESPythonGetType_PYTHON(PetscSNES snes, const char *name[]) \ 2097 except PETSC_ERR_PYTHON with gil: 2098 FunctionBegin(b"SNESPythonGetType_PYTHON") 2099 name[0] = PySNES(snes).getname() 2100 return FunctionEnd() 2101 2102cdef PetscErrorCode SNESCreate_Python( 2103 PetscSNES snes, 2104 ) except PETSC_ERR_PYTHON with gil: 2105 FunctionBegin(b"SNESCreate_Python") 2106 cdef SNESOps ops = snes.ops 2107 cdef PetscSNESLineSearch ls = NULL 2108 ops.reset = SNESReset_Python 2109 ops.destroy = SNESDestroy_Python 2110 ops.setup = SNESSetUp_Python 2111 ops.setfromoptions = SNESSetFromOptions_Python 2112 ops.view = SNESView_Python 2113 ops.solve = SNESSolve_Python 2114 # 2115 CHKERR(SNESParametersInitialize(snes)) 2116 CHKERR(PetscObjectComposeFunction( 2117 <PetscObject>snes, b"SNESPythonSetType_C", 2118 <PetscVoidFunction>SNESPythonSetType_PYTHON)) 2119 CHKERR(PetscObjectComposeFunction( 2120 <PetscObject>snes, b"SNESPythonGetType_C", 2121 <PetscVoidFunction>SNESPythonGetType_PYTHON)) 2122 # 2123 cdef ctx = PySNES(NULL) 2124 snes.data = <void*> ctx 2125 2126 # Ensure that the SNES has a linesearch object early enough that 2127 # it gets setFromOptions. 2128 CHKERR(SNESGetLineSearch(snes, &ls)) 2129 Py_INCREF(<PyObject*>snes.data) 2130 return FunctionEnd() 2131 2132cdef inline PetscErrorCode SNESDestroy_Python_inner( 2133 PetscSNES snes, 2134 ) except PETSC_ERR_PYTHON with gil: 2135 try: 2136 addRef(snes) 2137 SNESPythonSetContext(snes, NULL) 2138 finally: 2139 delRef(snes) 2140 Py_DECREF(<PyObject*>snes.data) 2141 snes.data = NULL 2142 return PETSC_SUCCESS 2143 2144cdef PetscErrorCode SNESDestroy_Python( 2145 PetscSNES snes, 2146 ) except PETSC_ERR_PYTHON nogil: 2147 FunctionBegin(b"SNESDestroy_Python") 2148 CHKERR(PetscObjectComposeFunction( 2149 <PetscObject>snes, b"SNESPythonSetType_C", 2150 <PetscVoidFunction>NULL)) 2151 CHKERR(PetscObjectComposeFunction( 2152 <PetscObject>snes, b"SNESPythonGetType_C", 2153 <PetscVoidFunction>NULL)) 2154 # 2155 if Py_IsInitialized(): SNESDestroy_Python_inner(snes) 2156 return FunctionEnd() 2157 2158cdef PetscErrorCode SNESSetUp_Python( 2159 PetscSNES snes, 2160 ) except PETSC_ERR_PYTHON with gil: 2161 FunctionBegin(b"SNESSetUp_Python") 2162 cdef char name[2048] 2163 cdef PetscBool found = PETSC_FALSE 2164 if PySNES(snes).self is None: 2165 CHKERR(PetscOptionsGetString(NULL, 2166 getPrefix(snes), b"-snes_python_type", 2167 name, sizeof(name), &found)) 2168 if found and name[0]: 2169 CHKERR(SNESPythonSetType_PYTHON(snes, name)) 2170 if PySNES(snes).self is None: 2171 return PetscSETERR(PETSC_ERR_USER, 2172 "Python context not set, call one of \n" 2173 " * SNESPythonSetType(snes, \"[package.]module.class\")\n" 2174 " * SNESSetFromOptions(snes) and pass option " 2175 "-snes_python_type [package.]module.class") 2176 # 2177 cdef setUp = PySNES(snes).setUp 2178 if setUp is not None: 2179 setUp(SNES_(snes)) 2180 return FunctionEnd() 2181 2182cdef inline PetscErrorCode SNESReset_Python_inner( 2183 PetscSNES snes, 2184 ) except PETSC_ERR_PYTHON with gil: 2185 cdef reset = PySNES(snes).reset 2186 if reset is not None: 2187 reset(SNES_(snes)) 2188 return PETSC_SUCCESS 2189 2190cdef PetscErrorCode SNESReset_Python( 2191 PetscSNES snes, 2192 ) except PETSC_ERR_PYTHON nogil: 2193 if getRef(snes) == 0: return PETSC_SUCCESS 2194 FunctionBegin(b"SNESReset_Python") 2195 if Py_IsInitialized(): SNESReset_Python_inner(snes) 2196 return FunctionEnd() 2197 2198cdef PetscErrorCode SNESSetFromOptions_Python( 2199 PetscSNES snes, 2200 PetscOptionItems PetscOptionsObject, 2201 ) except PETSC_ERR_PYTHON with gil: 2202 FunctionBegin(b"SNESSetFromOptions_Python") 2203 cdef char name[2048], *defval = PySNES(snes).getname() 2204 cdef PetscBool found = PETSC_FALSE 2205 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 2206 CHKERR(PetscOptionsString( 2207 b"-snes_python_type", b"Python [package.]module[.{class|function}]", 2208 b"SNESPythonSetType", defval, name, sizeof(name), &found)); <void>opts 2209 if found and name[0]: 2210 CHKERR(SNESPythonSetType_PYTHON(snes, name)) 2211 # 2212 cdef setFromOptions = PySNES(snes).setFromOptions 2213 if setFromOptions is not None: 2214 setFromOptions(SNES_(snes)) 2215 return FunctionEnd() 2216 2217cdef PetscErrorCode SNESView_Python( 2218 PetscSNES snes, 2219 PetscViewer vwr, 2220 ) except PETSC_ERR_PYTHON with gil: 2221 FunctionBegin(b"SNESView_Python") 2222 viewcontext(PySNES(snes), vwr) 2223 cdef view = PySNES(snes).view 2224 if view is not None: 2225 view(SNES_(snes), Viewer_(vwr)) 2226 return FunctionEnd() 2227 2228cdef PetscErrorCode SNESSolve_Python( 2229 PetscSNES snes, 2230 ) except PETSC_ERR_PYTHON with gil: 2231 FunctionBegin(b"SNESSolve_Python") 2232 cdef PetscVec b = NULL, x = NULL 2233 CHKERR(SNESGetRhs(snes, &b)) 2234 CHKERR(SNESGetSolution(snes, &x)) 2235 # 2236 snes.iter = 0 2237 # 2238 cdef solve = PySNES(snes).solve 2239 if solve is not None: 2240 solve(SNES_(snes), Vec_(b) if b != NULL else None, Vec_(x)) 2241 else: 2242 SNESSolve_Python_default(snes) 2243 # 2244 return FunctionEnd() 2245 2246cdef PetscErrorCode SNESSolve_Python_default( 2247 PetscSNES snes, 2248 ) except PETSC_ERR_PYTHON with gil: 2249 FunctionBegin(b"SNESSolve_Python_default") 2250 cdef PetscVec X=NULL, F=NULL, Y=NULL 2251 cdef PetscSNESLineSearch ls=NULL 2252 CHKERR(SNESGetSolution(snes, &X)) 2253 CHKERR(SNESGetFunction(snes, &F, NULL, NULL)) 2254 CHKERR(SNESGetSolutionUpdate(snes, &Y)) 2255 CHKERR(SNESGetLineSearch(snes, &ls)) 2256 # 2257 CHKERR(VecSet(Y, 0.0)) 2258 snes.ynorm = 0.0 2259 CHKERR(SNESComputeFunction(snes, X, F)) 2260 CHKERR(VecNorm(X, PETSC_NORM_2, &snes.xnorm)) 2261 CHKERR(VecNorm(F, PETSC_NORM_2, &snes.norm)) 2262 # 2263 cdef PetscInt lits=0 2264 CHKERR(SNESLogHistory(snes, snes.norm, lits)) 2265 CHKERR(SNESConverged(snes, snes.iter, snes.xnorm, snes.ynorm, snes.norm)) 2266 CHKERR(SNESMonitor(snes, snes.iter, snes.norm)) 2267 if snes.reason: 2268 return FunctionEnd() 2269 2270 for its from 0 <= its < snes.max_its: 2271 <void> its # unused 2272 SNESComputeUpdate(snes) 2273 SNESPreStep_Python(snes) 2274 # 2275 lits = -snes.linear_its 2276 SNESStep_Python(snes, X, F, Y) 2277 lits += snes.linear_its 2278 # 2279 CHKERR(SNESLineSearchApply(ls, X, F, NULL, Y)) 2280 CHKERR(SNESLineSearchGetNorms(ls, &snes.xnorm, &snes.norm, &snes.ynorm)) 2281 snes.iter += 1 2282 # 2283 SNESPostStep_Python(snes) 2284 CHKERR(SNESLogHistory(snes, snes.norm, lits)) 2285 CHKERR(SNESConverged(snes, snes.iter, snes.xnorm, snes.ynorm, snes.norm)) 2286 CHKERR(SNESMonitor(snes, snes.iter, snes.norm)) 2287 if snes.reason: break 2288 # 2289 return FunctionEnd() 2290 2291cdef PetscErrorCode SNESPreStep_Python( 2292 PetscSNES snes, 2293 ) except PETSC_ERR_PYTHON with gil: 2294 FunctionBegin(b"SNESPreStep_Python") 2295 cdef preStep = PySNES(snes).preStep 2296 if preStep is not None: 2297 preStep(SNES_(snes)) 2298 return FunctionEnd() 2299 2300cdef PetscErrorCode SNESPostStep_Python( 2301 PetscSNES snes, 2302 ) except PETSC_ERR_PYTHON with gil: 2303 FunctionBegin(b"SNESPostStep_Python") 2304 cdef postStep = PySNES(snes).postStep 2305 if postStep is not None: 2306 postStep(SNES_(snes)) 2307 return FunctionEnd() 2308 2309cdef PetscErrorCode SNESStep_Python( 2310 PetscSNES snes, 2311 PetscVec X, 2312 PetscVec F, 2313 PetscVec Y, 2314 ) except PETSC_ERR_PYTHON with gil: 2315 FunctionBegin(b"SNESStep_Python") 2316 cdef step = PySNES(snes).step 2317 if step is not None: 2318 step(SNES_(snes), Vec_(X), Vec_(F), Vec_(Y)) 2319 else: 2320 SNESStep_Python_default(snes, X, F, Y) 2321 return FunctionEnd() 2322 2323cdef PetscErrorCode SNESStep_Python_default( 2324 PetscSNES snes, 2325 PetscVec X, 2326 PetscVec F, 2327 PetscVec Y, 2328 ) except PETSC_ERR_PYTHON with gil: 2329 FunctionBegin(b"SNESStep_Python_default") 2330 cdef PetscMat J = NULL, P = NULL 2331 cdef PetscInt lits = 0 2332 CHKERR(SNESGetJacobian(snes, &J, &P, NULL, NULL)) 2333 CHKERR(SNESComputeJacobian(snes, X, J, P)) 2334 CHKERR(KSPSetOperators(snes.ksp, J, P)) 2335 CHKERR(KSPSolve(snes.ksp, F, Y)) 2336 CHKERR(KSPGetIterationNumber(snes.ksp, &lits)) 2337 snes.linear_its += lits 2338 return FunctionEnd() 2339 2340# -------------------------------------------------------------------- 2341 2342 2343cdef extern from * nogil: 2344 struct _TSOps: 2345 PetscErrorCode (*destroy)(PetscTS) except PETSC_ERR_PYTHON 2346 PetscErrorCode (*setup)(PetscTS) except PETSC_ERR_PYTHON 2347 PetscErrorCode (*reset)(PetscTS) except PETSC_ERR_PYTHON 2348 PetscErrorCode (*setfromoptions)(PetscTS, PetscOptionItems) except PETSC_ERR_PYTHON 2349 PetscErrorCode (*view)(PetscTS, PetscViewer) except PETSC_ERR_PYTHON 2350 PetscErrorCode (*step)(PetscTS) except PETSC_ERR_PYTHON 2351 PetscErrorCode (*rollback)(PetscTS) except PETSC_ERR_PYTHON 2352 PetscErrorCode (*interpolate)(PetscTS, PetscReal, PetscVec) except PETSC_ERR_PYTHON 2353 PetscErrorCode (*evaluatestep)(PetscTS, PetscInt, PetscVec, PetscBool*) except PETSC_ERR_PYTHON 2354 PetscErrorCode (*solve)(PetscTS) except PETSC_ERR_PYTHON 2355 PetscErrorCode (*snesfunction)(PetscSNES, PetscVec, PetscVec, PetscTS) except PETSC_ERR_PYTHON 2356 PetscErrorCode (*snesjacobian)(PetscSNES, PetscVec, PetscMat, PetscMat, PetscTS) except PETSC_ERR_PYTHON 2357 ctypedef _TSOps *TSOps 2358 struct _TSUserOps: 2359 PetscErrorCode (*prestep)(PetscTS) except PETSC_ERR_PYTHON 2360 PetscErrorCode (*prestage)(PetscTS, PetscReal) except PETSC_ERR_PYTHON 2361 PetscErrorCode (*poststage)(PetscTS, PetscReal, PetscInt, PetscVec*) except PETSC_ERR_PYTHON 2362 PetscErrorCode (*poststep)(PetscTS) except PETSC_ERR_PYTHON 2363 PetscErrorCode (*rhsfunction)(PetscTS, PetscReal, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 2364 PetscErrorCode (*ifunction) (PetscTS, PetscReal, PetscVec, PetscVec, PetscVec, void*) except PETSC_ERR_PYTHON 2365 PetscErrorCode (*rhsjacobian)(PetscTS, PetscReal, PetscVec, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 2366 PetscErrorCode (*ijacobian) (PetscTS, PetscReal, PetscVec, PetscVec, PetscReal, PetscMat, PetscMat, void*) except PETSC_ERR_PYTHON 2367 ctypedef _TSUserOps *TSUserOps 2368 struct _p_TS: 2369 void *data 2370 PetscDM dm 2371 PetscTSAdapt adapt 2372 TSOps ops 2373 TSUserOps userops 2374 PetscTSProblemType problem_type 2375 PetscInt snes_its 2376 PetscInt ksp_its 2377 PetscInt reject 2378 PetscInt max_reject 2379 PetscInt steps 2380 PetscReal ptime 2381 PetscVec vec_sol 2382 PetscReal time_step 2383 PetscInt max_steps 2384 PetscReal max_time 2385 PetscTSConvergedReason reason 2386 PetscSNES snes 2387 PetscBool usessnes 2388 2389 2390@cython.internal 2391cdef class _PyTS(_PyObj): pass 2392cdef inline _PyTS PyTS(PetscTS ts): 2393 if ts != NULL and ts.data != NULL: 2394 return <_PyTS>ts.data 2395 else: 2396 return _PyTS.__new__(_PyTS) 2397 2398cdef public PetscErrorCode TSPythonGetContext(PetscTS ts, void **ctx) \ 2399 except PETSC_ERR_PYTHON: 2400 FunctionBegin(b"TSPythonGetContext") 2401 PyTS(ts).getcontext(ctx) 2402 return FunctionEnd() 2403 2404cdef public PetscErrorCode TSPythonSetContext(PetscTS ts, void *ctx) \ 2405 except PETSC_ERR_PYTHON: 2406 FunctionBegin(b"TSPythonSetContext") 2407 PyTS(ts).setcontext(ctx, TS_(ts)) 2408 return FunctionEnd() 2409 2410cdef PetscErrorCode TSPythonSetType_PYTHON(PetscTS ts, const char *name) \ 2411 except PETSC_ERR_PYTHON with gil: 2412 FunctionBegin(b"TSPythonSetType_PYTHON") 2413 if name == NULL: return FunctionEnd() # XXX 2414 cdef object ctx = createcontext(name) 2415 TSPythonSetContext(ts, <void*>ctx) 2416 PyTS(ts).setname(name) 2417 return FunctionEnd() 2418 2419cdef PetscErrorCode TSPythonGetType_PYTHON(PetscTS ts, const char *name[]) \ 2420 except PETSC_ERR_PYTHON with gil: 2421 FunctionBegin(b"TSPythonGetType_PYTHON") 2422 name[0] = PyTS(ts).getname() 2423 return FunctionEnd() 2424 2425cdef PetscErrorCode TSCreate_Python( 2426 PetscTS ts, 2427 ) except PETSC_ERR_PYTHON with gil: 2428 FunctionBegin(b"TSCreate_Python") 2429 cdef TSOps ops = ts.ops 2430 ops.reset = TSReset_Python 2431 ops.destroy = TSDestroy_Python 2432 ops.setup = TSSetUp_Python 2433 ops.setfromoptions = TSSetFromOptions_Python 2434 ops.view = TSView_Python 2435 ops.step = TSStep_Python 2436 ops.rollback = TSRollBack_Python 2437 ops.interpolate = TSInterpolate_Python 2438 ops.evaluatestep = TSEvaluateStep_Python 2439 ops.snesfunction = SNESTSFormFunction_Python 2440 ops.snesjacobian = SNESTSFormJacobian_Python 2441 # 2442 CHKERR(PetscObjectComposeFunction( 2443 <PetscObject>ts, b"TSPythonSetType_C", 2444 <PetscVoidFunction>TSPythonSetType_PYTHON)) 2445 CHKERR(PetscObjectComposeFunction( 2446 <PetscObject>ts, b"TSPythonGetType_C", 2447 <PetscVoidFunction>TSPythonGetType_PYTHON)) 2448 # 2449 ts.usessnes = PETSC_TRUE 2450 # 2451 cdef ctx = PyTS(NULL) 2452 ts.data = <void*> ctx 2453 Py_INCREF(<PyObject*>ts.data) 2454 return FunctionEnd() 2455 2456cdef inline PetscErrorCode TSDestroy_Python_inner( 2457 PetscTS ts, 2458 ) except PETSC_ERR_PYTHON with gil: 2459 try: 2460 addRef(ts) 2461 TSPythonSetContext(ts, NULL) 2462 finally: 2463 delRef(ts) 2464 Py_DECREF(<PyObject*>ts.data) 2465 ts.data = NULL 2466 return PETSC_SUCCESS 2467 2468cdef PetscErrorCode TSDestroy_Python( 2469 PetscTS ts, 2470 ) except PETSC_ERR_PYTHON nogil: 2471 FunctionBegin(b"TSDestroy_Python") 2472 CHKERR(PetscObjectComposeFunction( 2473 <PetscObject>ts, b"TSPythonSetType_C", 2474 <PetscVoidFunction>NULL)) 2475 CHKERR(PetscObjectComposeFunction( 2476 <PetscObject>ts, b"TSPythonGetType_C", 2477 <PetscVoidFunction>NULL)) 2478 # 2479 if Py_IsInitialized(): TSDestroy_Python_inner(ts) 2480 return FunctionEnd() 2481 2482cdef PetscErrorCode TSSetUp_Python( 2483 PetscTS ts, 2484 ) except PETSC_ERR_PYTHON with gil: 2485 FunctionBegin(b"TSSetUp_Python") 2486 cdef PetscVec vec_update = NULL 2487 CHKERR(VecDuplicate(ts.vec_sol, &vec_update)) 2488 CHKERR(PetscObjectCompose(<PetscObject>ts, 2489 b"@ts.vec_update", 2490 <PetscObject>vec_update)) 2491 CHKERR(VecDestroy(&vec_update)) 2492 cdef PetscVec vec_dot = NULL 2493 CHKERR(VecDuplicate(ts.vec_sol, &vec_dot)) 2494 CHKERR(PetscObjectCompose(<PetscObject>ts, 2495 b"@ts.vec_dot", 2496 <PetscObject>vec_dot)) 2497 CHKERR(VecDestroy(&vec_dot)) 2498 # 2499 cdef char name[2048] 2500 cdef PetscBool found = PETSC_FALSE 2501 if PyTS(ts).self is None: 2502 CHKERR(PetscOptionsGetString(NULL, 2503 getPrefix(ts), b"-ts_python_type", 2504 name, sizeof(name), &found)) 2505 if found and name[0]: 2506 CHKERR(TSPythonSetType_PYTHON(ts, name)) 2507 if PyTS(ts).self is None: 2508 return PetscSETERR(PETSC_ERR_USER, 2509 "Python context not set, call one of \n" 2510 " * TSPythonSetType(ts, \"[package.]module.class\")\n" 2511 " * TSSetFromOptions(ts) and pass option " 2512 "-ts_python_type [package.]module.class") 2513 # 2514 cdef setUp = PyTS(ts).setUp 2515 if setUp is not None: 2516 setUp(TS_(ts)) 2517 return FunctionEnd() 2518 2519cdef inline PetscErrorCode TSReset_Python_inner( 2520 PetscTS ts, 2521 ) except PETSC_ERR_PYTHON with gil: 2522 cdef reset = PyTS(ts).reset 2523 if reset is not None: 2524 reset(TS_(ts)) 2525 return PETSC_SUCCESS 2526 2527cdef PetscErrorCode TSReset_Python( 2528 PetscTS ts, 2529 ) except PETSC_ERR_PYTHON nogil: 2530 if getRef(ts) == 0: return PETSC_SUCCESS 2531 FunctionBegin(b"TSReset_Python") 2532 CHKERR(PetscObjectCompose(<PetscObject>ts, b"@ts.vec_update", NULL)) 2533 CHKERR(PetscObjectCompose(<PetscObject>ts, b"@ts.vec_dot", NULL)) 2534 if Py_IsInitialized(): TSReset_Python_inner(ts) 2535 return FunctionEnd() 2536 2537cdef PetscErrorCode TSSetFromOptions_Python( 2538 PetscTS ts, 2539 PetscOptionItems PetscOptionsObject, 2540 ) except PETSC_ERR_PYTHON with gil: 2541 FunctionBegin(b"TSSetFromOptions_Python") 2542 cdef char name[2048], *defval = PyTS(ts).getname() 2543 cdef PetscBool found = PETSC_FALSE 2544 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 2545 CHKERR(PetscOptionsString( 2546 b"-ts_python_type", b"Python [package.]module[.{class|function}]", 2547 b"TSPythonSetType", defval, name, sizeof(name), &found)); <void>opts 2548 if found and name[0]: 2549 CHKERR(TSPythonSetType_PYTHON(ts, name)) 2550 # 2551 cdef setFromOptions = PyTS(ts).setFromOptions 2552 if setFromOptions is not None: 2553 setFromOptions(TS_(ts)) 2554 return FunctionEnd() 2555 2556cdef PetscErrorCode TSView_Python( 2557 PetscTS ts, 2558 PetscViewer vwr, 2559 ) except PETSC_ERR_PYTHON with gil: 2560 FunctionBegin(b"TSView_Python") 2561 viewcontext(PyTS(ts), vwr) 2562 cdef view = PyTS(ts).view 2563 if view is not None: 2564 view(TS_(ts), Viewer_(vwr)) 2565 return FunctionEnd() 2566 2567cdef PetscErrorCode TSStep_Python( 2568 PetscTS ts, 2569 ) except PETSC_ERR_PYTHON with gil: 2570 FunctionBegin(b"TSStep_Python") 2571 cdef step = PyTS(ts).step 2572 if step is not None: 2573 step(TS_(ts)) 2574 else: 2575 TSStep_Python_default(ts) 2576 return FunctionEnd() 2577 2578cdef PetscErrorCode TSRollBack_Python( 2579 PetscTS ts, 2580 ) except PETSC_ERR_PYTHON with gil: 2581 FunctionBegin(b"TSRollBack_Python") 2582 cdef rollback = PyTS(ts).rollback 2583 if rollback is None: return UNSUPPORTED(b"rollback") 2584 rollback(TS_(ts)) 2585 return FunctionEnd() 2586 2587cdef PetscErrorCode TSInterpolate_Python( 2588 PetscTS ts, 2589 PetscReal t, 2590 PetscVec x, 2591 ) except PETSC_ERR_PYTHON with gil: 2592 FunctionBegin(b"TSInterpolate _Python") 2593 cdef interpolate = PyTS(ts).interpolate 2594 if interpolate is None: return UNSUPPORTED(b"interpolate") 2595 interpolate(TS_(ts), toReal(t), Vec_(x)) 2596 return FunctionEnd() 2597 2598cdef PetscErrorCode TSEvaluateStep_Python( 2599 PetscTS ts, 2600 PetscInt o, 2601 PetscVec x, 2602 PetscBool *flag, 2603 ) except PETSC_ERR_PYTHON with gil: 2604 FunctionBegin(b"TSEvaluateStep _Python") 2605 cdef evaluatestep = PyTS(ts).evaluatestep 2606 if evaluatestep is None: return UNSUPPORTED(b"evaluatestep") 2607 cdef done = evaluatestep(TS_(ts), toInt(o), Vec_(x)) 2608 if flag != NULL: 2609 flag[0] = PETSC_TRUE if done else PETSC_FALSE 2610 elif not done: 2611 return PetscSETERR(PETSC_ERR_USER, "Cannot evaluate step") 2612 return FunctionEnd() 2613 2614cdef PetscErrorCode SNESTSFormFunction_Python( 2615 PetscSNES snes, 2616 PetscVec x, 2617 PetscVec f, 2618 PetscTS ts, 2619 ) except PETSC_ERR_PYTHON with gil: 2620 FunctionBegin(b"SNESTSFormFunction _Python") 2621 cdef formSNESFunction = PyTS(ts).formSNESFunction 2622 if formSNESFunction is not None: 2623 args = (SNES_(snes), Vec_(x), Vec_(f), TS_(ts)) 2624 formSNESFunction(args) 2625 return FunctionEnd() 2626 # 2627 cdef PetscVec dx = NULL 2628 CHKERR(PetscObjectQuery( 2629 <PetscObject>ts, 2630 b"@ts.vec_dot", 2631 <PetscObject*>&dx)) 2632 # 2633 cdef PetscReal t = ts.ptime + ts.time_step 2634 cdef PetscReal a = 1.0/ts.time_step 2635 CHKERR(VecCopy(ts.vec_sol, dx)) 2636 CHKERR(VecAXPBY(dx, +a, -a, x)) 2637 CHKERR(TSComputeIFunction(ts, t, x, dx, f, PETSC_FALSE)) 2638 return FunctionEnd() 2639 2640cdef PetscErrorCode SNESTSFormJacobian_Python( 2641 PetscSNES snes, 2642 PetscVec x, 2643 PetscMat A, 2644 PetscMat B, 2645 PetscTS ts, 2646 ) except PETSC_ERR_PYTHON with gil: 2647 FunctionBegin(b"SNESTSFormJacobian _Python") 2648 cdef formSNESJacobian = PyTS(ts).formSNESJacobian 2649 if formSNESJacobian is not None: 2650 args = (SNES_(snes), Vec_(x), Mat_(A), Mat_(B), TS_(ts)) 2651 formSNESJacobian(*args) 2652 return FunctionEnd() 2653 # 2654 cdef PetscVec dx = NULL 2655 CHKERR(PetscObjectQuery( 2656 <PetscObject>ts, 2657 b"@ts.vec_dot", 2658 <PetscObject*>&dx)) 2659 # 2660 cdef PetscReal t = ts.ptime + ts.time_step 2661 cdef PetscReal a = 1.0/ts.time_step 2662 CHKERR(VecCopy(ts.vec_sol, dx)) 2663 CHKERR(VecAXPBY(dx, +a, -a, x)) 2664 CHKERR(TSComputeIJacobian(ts, t, x, dx, a, A, B, PETSC_FALSE)) 2665 return FunctionEnd() 2666 2667cdef PetscErrorCode TSSolveStep_Python( 2668 PetscTS ts, 2669 PetscReal t, 2670 PetscVec x, 2671 ) except PETSC_ERR_PYTHON with gil: 2672 FunctionBegin(b"TSSolveStep_Python") 2673 cdef solveStep = PyTS(ts).solveStep 2674 if solveStep is not None: 2675 solveStep(TS_(ts), <double>t, Vec_(x)) 2676 return FunctionEnd() 2677 # 2678 cdef PetscInt nits = 0, lits = 0 2679 CHKERR(SNESSolve(ts.snes, NULL, x)) 2680 CHKERR(SNESGetIterationNumber(ts.snes, &nits)) 2681 CHKERR(SNESGetLinearSolveIterations(ts.snes, &lits)) 2682 ts.snes_its += nits 2683 ts.ksp_its += lits 2684 return FunctionEnd() 2685 2686cdef PetscErrorCode TSAdaptStep_Python( 2687 PetscTS ts, 2688 PetscReal t, 2689 PetscVec x, 2690 PetscReal *nextdt, 2691 PetscBool *stepok, 2692 ) except PETSC_ERR_PYTHON with gil: 2693 FunctionBegin(b"TSAdaptStep_Python") 2694 nextdt[0] = ts.time_step 2695 stepok[0] = PETSC_TRUE 2696 # 2697 cdef adaptStep = PyTS(ts).adaptStep 2698 if adaptStep is None: return FunctionEnd() 2699 cdef object retval 2700 cdef double dt 2701 cdef bint ok 2702 retval = adaptStep(TS_(ts), <double>t, Vec_(x)) 2703 if retval is None: pass 2704 elif isinstance(retval, float): 2705 dt = retval 2706 nextdt[0] = <PetscReal>dt 2707 stepok[0] = PETSC_TRUE 2708 elif isinstance(retval, bool): 2709 ok = retval 2710 nextdt[0] = ts.time_step 2711 stepok[0] = PETSC_TRUE if ok else PETSC_FALSE 2712 else: 2713 dt, ok = retval 2714 nextdt[0] = <PetscReal>dt 2715 stepok[0] = PETSC_TRUE if ok else PETSC_FALSE 2716 return FunctionEnd() 2717 2718cdef PetscErrorCode TSStep_Python_default( 2719 PetscTS ts, 2720 ) except PETSC_ERR_PYTHON with gil: 2721 FunctionBegin(b"TSStep_Python_default") 2722 cdef PetscVec vec_update = NULL 2723 CHKERR(PetscObjectQuery( 2724 <PetscObject>ts, 2725 b"@ts.vec_update", 2726 <PetscObject*>&vec_update)) 2727 # 2728 cdef PetscReal tt = ts.ptime 2729 cdef PetscReal dt = ts.time_step 2730 cdef PetscBool accept = PETSC_TRUE 2731 cdef PetscBool stageok = PETSC_TRUE 2732 for r from 0 <= r < ts.max_reject: 2733 <void> r # unused 2734 tt = ts.ptime + ts.time_step 2735 CHKERR(VecCopy(ts.vec_sol, vec_update)) 2736 CHKERR(TSPreStage(ts, tt+dt)) 2737 TSSolveStep_Python(ts, tt, vec_update) 2738 CHKERR(TSPostStage(ts, tt+dt, 0, &vec_update)) 2739 CHKERR(TSAdaptCheckStage(ts.adapt, ts, tt+dt, vec_update, &stageok)) 2740 if not stageok: 2741 ts.reject += 1 2742 continue 2743 TSAdaptStep_Python(ts, tt, vec_update, &dt, &accept) 2744 if not accept: 2745 ts.time_step = dt 2746 ts.reject += 1 2747 continue 2748 CHKERR(VecCopy(vec_update, ts.vec_sol)) 2749 ts.ptime += ts.time_step 2750 ts.time_step = dt 2751 break 2752 if (not stageok or not accept) and ts.reason == 0: 2753 ts.reason = TS_DIVERGED_STEP_REJECTED 2754 return FunctionEnd() 2755 2756# -------------------------------------------------------------------- 2757 2758 2759cdef extern from * nogil: 2760 struct _TaoOps: 2761 PetscErrorCode (*destroy)(PetscTAO) except PETSC_ERR_PYTHON 2762 PetscErrorCode (*setup)(PetscTAO) except PETSC_ERR_PYTHON 2763 PetscErrorCode (*solve)(PetscTAO) except PETSC_ERR_PYTHON 2764 PetscErrorCode (*setfromoptions)(PetscTAO, PetscOptionItems) except PETSC_ERR_PYTHON 2765 PetscErrorCode (*view)(PetscTAO, PetscViewer) except PETSC_ERR_PYTHON 2766 ctypedef _TaoOps *TaoOps 2767 struct _p_TAO: 2768 void *data 2769 TaoOps ops 2770 PetscInt niter, max_it 2771 PetscTAOConvergedReason reason 2772 PetscInt ksp_its, ksp_tot_its 2773 PetscKSP ksp 2774 PetscVec gradient 2775 PetscVec stepdirection 2776 PetscTAOLineSearch linesearch 2777 2778cdef extern from * nogil: # custom.h 2779 PetscErrorCode TaoConverged(PetscTAO, PetscTAOConvergedReason*) 2780 2781cdef extern from * nogil: # custom.h 2782 PetscErrorCode TaoGetVecs(PetscTAO, PetscVec*, PetscVec*, PetscVec*) 2783 PetscErrorCode TaoCheckReals(PetscTAO, PetscReal, PetscReal) 2784 PetscErrorCode TaoComputeUpdate(PetscTAO, PetscReal*) 2785 PetscErrorCode TaoCreateDefaultLineSearch(PetscTAO) 2786 PetscErrorCode TaoCreateDefaultKSP(PetscTAO) 2787 PetscErrorCode TaoApplyLineSearch(PetscTAO, PetscReal*, PetscReal*, PetscTAOLineSearchConvergedReason*) 2788 2789 2790@cython.internal 2791cdef class _PyTao(_PyObj): pass 2792cdef inline _PyTao PyTao(PetscTAO tao): 2793 if tao != NULL and tao.data != NULL: 2794 return <_PyTao>tao.data 2795 else: 2796 return _PyTao.__new__(_PyTao) 2797 2798cdef public PetscErrorCode TaoPythonGetContext(PetscTAO tao, void **ctx) \ 2799 except PETSC_ERR_PYTHON: 2800 FunctionBegin(b"TaoPythonGetContext") 2801 PyTao(tao).getcontext(ctx) 2802 return FunctionEnd() 2803 2804cdef public PetscErrorCode TaoPythonSetContext(PetscTAO tao, void *ctx) \ 2805 except PETSC_ERR_PYTHON: 2806 FunctionBegin(b"TaoPythonSetContext") 2807 PyTao(tao).setcontext(ctx, TAO_(tao)) 2808 return FunctionEnd() 2809 2810cdef PetscErrorCode TaoPythonSetType_PYTHON(PetscTAO tao, const char *name) \ 2811 except PETSC_ERR_PYTHON with gil: 2812 FunctionBegin(b"TaoPythonSetType_PYTHON") 2813 if name == NULL: return FunctionEnd() # XXX 2814 cdef object ctx = createcontext(name) 2815 TaoPythonSetContext(tao, <void*>ctx) 2816 PyTao(tao).setname(name) 2817 return FunctionEnd() 2818 2819cdef PetscErrorCode TaoPythonGetType_PYTHON(PetscTAO tao, const char *name[]) \ 2820 except PETSC_ERR_PYTHON with gil: 2821 FunctionBegin(b"TaoPythonGetType_PYTHON") 2822 name[0] = PyTao(tao).getname() 2823 return FunctionEnd() 2824 2825cdef PetscErrorCode TaoCreate_Python( 2826 PetscTAO tao, 2827 ) except PETSC_ERR_PYTHON with gil: 2828 FunctionBegin(b"TaoCreate_Python") 2829 cdef TaoOps ops = tao.ops 2830 ops.destroy = TaoDestroy_Python 2831 ops.view = TaoView_Python 2832 ops.solve = TaoSolve_Python 2833 ops.setup = TaoSetUp_Python 2834 ops.setfromoptions = TaoSetFromOptions_Python 2835 # 2836 CHKERR(TaoParametersInitialize(tao)) 2837 CHKERR(PetscObjectComposeFunction( 2838 <PetscObject>tao, b"TaoPythonSetType_C", 2839 <PetscVoidFunction>TaoPythonSetType_PYTHON)) 2840 CHKERR(PetscObjectComposeFunction( 2841 <PetscObject>tao, b"TaoPythonGetType_C", 2842 <PetscVoidFunction>TaoPythonGetType_PYTHON)) 2843 # 2844 CHKERR(TaoCreateDefaultLineSearch(tao)) 2845 CHKERR(TaoCreateDefaultKSP(tao)) 2846 # 2847 cdef ctx = PyTao(NULL) 2848 tao.data = <void*> ctx 2849 Py_INCREF(<PyObject*>tao.data) 2850 return FunctionEnd() 2851 2852cdef inline PetscErrorCode TaoDestroy_Python_inner( 2853 PetscTAO tao, 2854 ) except PETSC_ERR_PYTHON with gil: 2855 try: 2856 addRef(tao) 2857 TaoPythonSetContext(tao, NULL) 2858 finally: 2859 delRef(tao) 2860 Py_DECREF(<PyObject*>tao.data) 2861 tao.data = NULL 2862 return PETSC_SUCCESS 2863 2864cdef PetscErrorCode TaoDestroy_Python( 2865 PetscTAO tao, 2866 ) except PETSC_ERR_PYTHON nogil: 2867 FunctionBegin(b"TaoDestroy_Python") 2868 CHKERR(PetscObjectComposeFunction( 2869 <PetscObject>tao, b"TaoPythonSetType_C", 2870 <PetscVoidFunction>NULL)) 2871 CHKERR(PetscObjectComposeFunction( 2872 <PetscObject>tao, b"TaoPythonGetType_C", 2873 <PetscVoidFunction>NULL)) 2874 # 2875 if Py_IsInitialized(): TaoDestroy_Python_inner(tao) 2876 return FunctionEnd() 2877 2878cdef PetscErrorCode TaoSetUp_Python( 2879 PetscTAO tao, 2880 ) except PETSC_ERR_PYTHON with gil: 2881 FunctionBegin(b"TaoSetUp_Python") 2882 cdef char name[2048] 2883 cdef PetscBool found = PETSC_FALSE 2884 if PyTao(tao).self is None: 2885 CHKERR(PetscOptionsGetString(NULL, 2886 getPrefix(tao), b"-tao_python_type", 2887 name, sizeof(name), &found)) 2888 if found and name[0]: 2889 CHKERR(TaoPythonSetType_PYTHON(tao, name)) 2890 if PyTao(tao).self is None: 2891 return PetscSETERR(PETSC_ERR_USER, 2892 "Python context not set, call one of \n" 2893 " * TaoPythonSetType(tao, \"[package.]module.class\")\n" 2894 " * TaoSetFromOptions(tao) and pass option " 2895 "-tao_python_type [package.]module.class") 2896 # 2897 cdef setUp = PyTao(tao).setUp 2898 if setUp is not None: 2899 setUp(TAO_(tao)) 2900 return FunctionEnd() 2901 2902cdef PetscErrorCode TaoSetFromOptions_Python( 2903 PetscTAO tao, 2904 PetscOptionItems PetscOptionsObject, 2905 ) except PETSC_ERR_PYTHON with gil: 2906 FunctionBegin(b"TaoSetFromOptions_Python") 2907 cdef char name[2048], *defval = PyTao(tao).getname() 2908 cdef PetscBool found = PETSC_FALSE 2909 cdef PetscOptionItems opts "PetscOptionsObject" = PetscOptionsObject 2910 CHKERR(PetscOptionsString( 2911 b"-tao_python_type", b"Python [package.]module[.{class|function}]", 2912 b"TaoPythonSetType", defval, name, sizeof(name), &found)); <void>opts 2913 if found and name[0]: 2914 CHKERR(TaoPythonSetType_PYTHON(tao, name)) 2915 # 2916 cdef setFromOptions = PyTao(tao).setFromOptions 2917 if setFromOptions is not None: 2918 setFromOptions(TAO_(tao)) 2919 CHKERR(KSPSetFromOptions(tao.ksp)) 2920 return FunctionEnd() 2921 2922cdef PetscErrorCode TaoView_Python( 2923 PetscTAO tao, 2924 PetscViewer vwr, 2925 ) except PETSC_ERR_PYTHON with gil: 2926 FunctionBegin(b"TaoView_Python") 2927 viewcontext(PyTao(tao), vwr) 2928 cdef view = PyTao(tao).view 2929 if view is not None: 2930 view(TAO_(tao), Viewer_(vwr)) 2931 return FunctionEnd() 2932 2933cdef PetscErrorCode TaoSolve_Python( 2934 PetscTAO tao, 2935 ) except PETSC_ERR_PYTHON with gil: 2936 FunctionBegin(b"TaoSolve_Python") 2937 tao.niter = 0 2938 tao.ksp_its = 0 2939 tao.reason = TAO_CONTINUE_ITERATING 2940 # 2941 cdef solve = PyTao(tao).solve 2942 if solve is not None: 2943 solve(TAO_(tao)) 2944 else: 2945 TaoSolve_Python_default(tao) 2946 # 2947 return FunctionEnd() 2948 2949cdef PetscErrorCode TaoSolve_Python_default( 2950 PetscTAO tao, 2951 ) except PETSC_ERR_PYTHON with gil: 2952 FunctionBegin(b"TaoSolve_Python_default") 2953 cdef PetscVec X = NULL, G = NULL, S = NULL 2954 CHKERR(TaoGetVecs(tao, &X, &G, &S)) 2955 # 2956 cdef PetscReal f = 0.0 2957 cdef PetscReal gnorm = 0.0 2958 cdef PetscReal step = 1.0 2959 # 2960 if G != NULL: 2961 CHKERR(TaoComputeObjectiveAndGradient(tao, X, &f, G)) 2962 CHKERR(VecNorm(G, PETSC_NORM_2, &gnorm)) 2963 else: 2964 CHKERR(TaoComputeObjective(tao, X, &f)) 2965 CHKERR(TaoCheckReals(tao, f, gnorm)) 2966 2967 CHKERR(TaoLogConvergenceHistory(tao, f, gnorm, 0.0, tao.ksp_its)) 2968 CHKERR(TaoMonitor(tao, tao.niter, f, gnorm, 0.0, step)) 2969 CHKERR(TaoConverged(tao, &tao.reason)) 2970 2971 cdef PetscTAOLineSearchConvergedReason lsr = TAOLINESEARCH_SUCCESS 2972 for its from 0 <= its < tao.max_it: 2973 <void> its # unused 2974 if tao.reason: break 2975 CHKERR(TaoComputeUpdate(tao, &f)) 2976 2977 TaoPreStep_Python(tao) 2978 # 2979 tao.ksp_its = 0 2980 TaoStep_Python(tao, X, G, S) 2981 CHKERR(KSPGetIterationNumber(tao.ksp, &tao.ksp_its)) 2982 tao.ksp_tot_its += tao.ksp_its 2983 # 2984 if G != NULL: 2985 CHKERR(TaoApplyLineSearch(tao, &f, &step, &lsr)) 2986 CHKERR(VecNorm(G, PETSC_NORM_2, &gnorm)) 2987 if lsr < TAOLINESEARCH_CONTINUE_ITERATING: 2988 tao.reason = TAO_DIVERGED_LS_FAILURE 2989 else: 2990 CHKERR(TaoComputeObjective(tao, X, &f)) 2991 CHKERR(TaoCheckReals(tao, f, gnorm)) 2992 2993 tao.niter += 1 2994 # 2995 TaoPostStep_Python(tao) 2996 CHKERR(TaoLogConvergenceHistory(tao, f, gnorm, 0.0, tao.ksp_its)) 2997 CHKERR(TaoMonitor(tao, tao.niter, f, gnorm, 0.0, step)) 2998 CHKERR(TaoConverged(tao, &tao.reason)) 2999 3000 if tao.niter == tao.max_it: 3001 if tao.reason <= 0: 3002 tao.reason = TAO_DIVERGED_MAXITS 3003 # 3004 return FunctionEnd() 3005 3006cdef PetscErrorCode TaoStep_Python( 3007 PetscTAO tao, 3008 PetscVec X, 3009 PetscVec G, 3010 PetscVec S, 3011 ) except PETSC_ERR_PYTHON with gil: 3012 FunctionBegin(b"TaoStep_Python") 3013 cdef step = PyTao(tao).step 3014 if step is not None: 3015 step(TAO_(tao), Vec_(X), Vec_(G) if G != NULL else None, Vec_(S) if S != NULL else None) 3016 else: 3017 # TaoStep_Python_default(tao, X, G, S) 3018 CHKERR(TaoComputeGradient(tao, X, G)) 3019 CHKERR(VecAXPBY(S, -1.0, 0.0, G)) 3020 return FunctionEnd() 3021 3022cdef PetscErrorCode TaoPreStep_Python( 3023 PetscTAO tao, 3024 ) except PETSC_ERR_PYTHON with gil: 3025 FunctionBegin(b"TaoPreStep_Python") 3026 cdef preStep = PyTao(tao).preStep 3027 if preStep is not None: 3028 preStep(TAO_(tao)) 3029 return FunctionEnd() 3030 3031cdef PetscErrorCode TaoPostStep_Python( 3032 PetscTAO tao, 3033 ) except PETSC_ERR_PYTHON with gil: 3034 FunctionBegin(b"TaoPostStep_Python") 3035 cdef postStep = PyTao(tao).postStep 3036 if postStep is not None: 3037 postStep(TAO_(tao)) 3038 return FunctionEnd() 3039 3040# -------------------------------------------------------------------- 3041 3042cdef PetscErrorCode PetscPythonMonitorSet_Python( 3043 PetscObject obj_p, 3044 const char *url_p, 3045 ) except PETSC_ERR_PYTHON with gil: 3046 FunctionBegin(b"PetscPythonMonitorSet_Python") 3047 assert obj_p != NULL 3048 assert url_p != NULL 3049 assert url_p[0] != 0 3050 # 3051 cdef PetscClassId classid = 0 3052 CHKERR(PetscObjectGetClassId(obj_p, &classid)) 3053 cdef type klass = PyPetscType_Lookup(classid) 3054 cdef Object ob = klass() 3055 ob.obj[0] = newRef(obj_p) 3056 # 3057 cdef url = bytes2str(url_p) 3058 if ':' in url: 3059 path, names = parse_url(url) 3060 else: 3061 path, names = url, 'monitor' 3062 module = load_module(path) 3063 for attr in names.split(', '): 3064 monitor = getattr(module, attr) 3065 if isinstance(monitor, type): 3066 monitor = monitor(ob) 3067 ob.setMonitor(monitor) 3068 # 3069 return FunctionEnd() 3070 3071# -------------------------------------------------------------------- 3072 3073cdef extern from * nogil: 3074 3075 ctypedef PetscErrorCode MatCreateFunction (PetscMat) except PETSC_ERR_PYTHON 3076 ctypedef PetscErrorCode PCCreateFunction (PetscPC) except PETSC_ERR_PYTHON 3077 ctypedef PetscErrorCode KSPCreateFunction (PetscKSP) except PETSC_ERR_PYTHON 3078 ctypedef PetscErrorCode SNESCreateFunction (PetscSNES) except PETSC_ERR_PYTHON 3079 ctypedef PetscErrorCode TSCreateFunction (PetscTS) except PETSC_ERR_PYTHON 3080 ctypedef PetscErrorCode TaoCreateFunction (PetscTAO) except PETSC_ERR_PYTHON 3081 ctypedef PetscErrorCode PetscViewerCreateFunction (PetscViewer) except PETSC_ERR_PYTHON 3082 3083 PetscErrorCode MatRegister (const char[], MatCreateFunction*) 3084 PetscErrorCode PCRegister (const char[], PCCreateFunction*) 3085 PetscErrorCode KSPRegister (const char[], KSPCreateFunction*) 3086 PetscErrorCode SNESRegister (const char[], SNESCreateFunction*) 3087 PetscErrorCode TSRegister (const char[], TSCreateFunction*) 3088 PetscErrorCode TaoRegister (const char[], TaoCreateFunction*) 3089 PetscErrorCode PetscViewerRegister (const char[], PetscViewerCreateFunction*) 3090 3091 PetscErrorCode (*PetscPythonMonitorSet_C) \ 3092 (PetscObject, const char[]) except PETSC_ERR_PYTHON 3093 3094 3095cdef public PetscErrorCode PetscPythonRegisterAll() except PETSC_ERR_PYTHON: 3096 FunctionBegin(b"PetscPythonRegisterAll") 3097 3098 # Python subtypes 3099 CHKERR(MatRegister(MATPYTHON, MatCreate_Python)) 3100 CHKERR(PCRegister(PCPYTHON, PCCreate_Python)) 3101 CHKERR(KSPRegister(KSPPYTHON, KSPCreate_Python)) 3102 CHKERR(SNESRegister(SNESPYTHON, SNESCreate_Python)) 3103 CHKERR(TSRegister(TSPYTHON, TSCreate_Python)) 3104 CHKERR(TaoRegister(TAOPYTHON, TaoCreate_Python)) 3105 CHKERR(PetscViewerRegister(PETSCVIEWERPYTHON, PetscViewerCreate_Python)) 3106 3107 # Python monitors 3108 global PetscPythonMonitorSet_C 3109 PetscPythonMonitorSet_C = PetscPythonMonitorSet_Python 3110 3111 return FunctionEnd() 3112 3113# -------------------------------------------------------------------- 3114