-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 On 18/02/15 12:34, Eike Mueller wrote:
Hi Lawrence,
thanks, updating PETSc, petsc4py, ufl and firedrake fixes this, i.e. I can project the expression. However, I still have an issue when using PETSc vectors twice.
What I do in my code is this:
(in the constructor) self._x = PETSc.Vec() self._x.create() self._x.setSizes((self._ndof, None)) self._x.setFromOptions() self._y = self._x.duplicate()
(in the solve routine) [line 208]: self._ksp.solve(self._y,self._x)
When I call the solve routine a second time (I do a warmup solve first) I still get the error that the vector is in the wrong state (see below). Do I have to unlock it?
Thanks,
Eike
File "driver.py", line 478, in <module> main(parameter_filename) File "driver.py", line 388, in main u,p,b = gravitywave_solver_matrixfree.solve(r_u,r_p,r_b) File "/Users/eikemueller/PostDocBath/EllipticSolvers/Firedrake_workspace/firedrake-helmholtzsolver/source/gravitywaves.py", line 208, in solve self._ksp.solve(self._y,self._x) File "KSP.pyx", line 353, in petsc4py.PETSc.KSP.solve (src/petsc4py.PETSc.c:139671) File "libpetsc4py.pyx", line 1330, in libpetsc4py.PCApply_Python (src/libpetsc4py/libpetsc4py.c:13951) File "/Users/eikemueller/PostDocBath/EllipticSolvers/Firedrake_workspace/firedrake-helmholtzsolver/source/mixedpreconditioners.py", line 138, in apply self._mixedarray.split(x,u,p) File "/Users/eikemueller/PostDocBath/EllipticSolvers/Firedrake_workspace/firedrake-helmholtzsolver/source/mixedarray.py", line 65, in split x.array[:] = v.array[min_idx:max_idx] File "Vec.pyx", line 815, in petsc4py.PETSc.Vec.array.__get__ (src/petsc4py.PETSc.c:91905) File "arraynpy.pxi", line 67, in petsc4py.PETSc.asarray (src/petsc4py.PETSc.c:7333) File "petscvec.pxi", line 469, in petsc4py.PETSc._Vec_buffer.__getreadbuffer__ (src/petsc4py.PETSc.c:19630) File "petscvec.pxi", line 459, in petsc4py.PETSc._Vec_buffer.getbuffer (src/petsc4py.PETSc.c:19460) petsc4py.PETSc.Error: error code 73 [0] KSPSolve() line 572 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/src/ksp/ksp/interface/itfunc.c
[0] KSPSolve_GMRES() line 234 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/src/ksp/ksp/impls/gmres/gmres.c
[0] KSPInitialResidual() line 63 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/src/ksp/ksp/interface/itres.c
[0] KSP_PCApply() line 235 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/include/petsc-private/kspimpl.h
[0] PCApply() line 449 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/src/ksp/pc/interface/precon.c
[0] VecGetArray() line 1646 in /Users/eikemueller/PostDocBath/EllipticSolvers/petsc/src/vec/vec/interface/rvector.c
[0] Object is in wrong state [0] Vec is locked read only, argument # 1
Aha, so the problem here is the line: x.array[:] = v.array[...] v is "locked", but petsc4py doesn't give you a way of doing VecGetArrayRead (which you're allowed to do in this circumstance), only VecGetArray (which you are not). Workaround: Make some ISes. PETSc.IS().createStride(max_idx - min_idx, start=global_min_idx, step=1, comm=PETSc.COMM_SELF) where global_min_idx is the global (not local) number of the minimum index on this process. You can compute it by doing: v.owner_range[0] + min_idx Now in the loop you do: for i, x in enumerate(args): iset = self.range(i) tmp = v.getSubVector(iset) v.copy(x) v.restoreSubVector(iset, tmp) Cheers, Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 iQEcBAEBAgAGBQJU5IrtAAoJECOc1kQ8PEYvDeoH/0ZNgQ+mFuGjBuj4jaJXoglQ M24W7pP1auwwrG1bHZvkmFYWogYptJKD3zvUzbl4tbEswkSxgRTKEsVz0cU9qP5h 3p5ezKKjszO43QlqY04mEpbQPRXTw0HPyudf4YB2T1FtFl/nUTaGk/DsBDvApTI8 qEUQ2+xALU1Ry2wIFGUedSNkz7mde/3k4WZo/PgGfJziPRbgguQWTr6ethaG6ywU nP1wA/9xGswiUK3z8IlWbg9eVIb997OGkBfw/6PeJpFXfAwbG2TfbSVTFZ8ACfV/ h2Y8XvhYWk9V1okYpmvEYNiRSAz4M86ekrjr38GLtcGT6LT2Hg/6pYE5mDozsiA= =A7ap -----END PGP SIGNATURE-----