Closed ctjacobs closed 9 years ago
Direct solvers are not supported on mixed systems by PETSc afair, so it's not a Firedrake limitation - correct me if that's wrong, @wence-.
It's not possible to use a direct solve with a MatNest, I have a PyOP2 branch that allows assembly into a monolithic matrix which would allow a direct solve (baij-matrices). Unfortunately it doesn't work with vector function spaces due to a petsc bug that I have reported without any response.
Is the patch for this in petsc and petsc4py MAPDES-master? I just tried giving this a go with up-to-date packages, and locally merged master
into the baij-matrices
PyOP2 branch. I got a load of conflicts which I tried to fix manually. I'm still getting an error in petsc4py though:
Traceback (most recent call last):
File "test.py", line 51, in <module>
solver = NonlinearVariationalSolver(problem, solver_parameters={'ksp_type':'preonly', 'pc_type':'lu'})
File "/data/firedrake/firedrake/variational_solver.py", line 115, in __init__
ctx.set_jacobian(self.snes)
File "/data/firedrake/firedrake/solving_utils.py", line 132, in set_jacobian
snes.setJacobian(self.form_jacobian, J=self._jacs[-1]._M.handle,
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 753, in handle
self._init()
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 351, in _init
self._init_nest()
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 548, in _init_nest
mat.createNest([[m.handle for m in row_] for row_ in self._blocks],
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 753, in handle
self._init()
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 355, in _init
self._init_block()
File "/usr/local/lib/python2.7/dist-packages/PyOP2-0.11.0_179_ga0030bb_dirty-py2.7-linux-x86_64.egg/pyop2/petsc_base.py", line 566, in _init_block
(self.sparsity._rowptr, self.sparsity._colidx, self._array))
File "PETSc/Mat.pyx", line 343, in petsc4py.PETSc.Mat.createAIJWithArrays (src/petsc4py.PETSc.c:109338)
ValueError: A matrix with 50 rows requires a row pointer of length 51 (given: 26)
There is a dangling change in petsc (not merged) and two open changes in pyop2 and firedrake, either of those give details.
Should now work, pass nest=False to the NLVProblem constructor (or nest=False to a solve(...) call) or, if you want to set the default parameters["matnest"] = False
.
Works for me. Thanks.
I'm trying to use a direct solver with a mixed problem and the solution to this issue no longer seems correct.
firedrake.LinearVariationalProblem.__init__
does not have a nest
keyword argument. Neither does firedrake.LinearVariationalSolver.solve
. The same is true for the nonlinear counterparts. Setting firedrake.parameters["matnest"] = False
does not change the behavior.
Here's a MWE based on the mixed Poisson tutorial:
import firedrake as fe
print("Using firedrake-" + fe.__version__)
mesh = fe.UnitSquareMesh(3, 3)
BDM = fe.FunctionSpace(mesh, "BDM", 1)
DG = fe.FunctionSpace(mesh, "DG", 0)
W = BDM * DG
sigma, u = fe.TrialFunctions(W)
tau, v = fe.TestFunctions(W)
x, y = fe.SpatialCoordinate(mesh)
f = fe.Function(DG).interpolate(
10*fe.exp(-(pow(x - 0.5, 2) + pow(y - 0.5, 2)) / 0.02))
dot, div, dx = fe.dot, fe.div, fe.dx
a = (dot(sigma, tau) + div(tau)*u + div(sigma)*v)*dx
L = - f*v*dx
bc0 = fe.DirichletBC(W.sub(0), fe.as_vector([0.0, -fe.sin(5*x)]), 1)
bc1 = fe.DirichletBC(W.sub(0), fe.as_vector([0.0, fe.sin(5*y)]), 2)
w = fe.Function(W)
fe.parameters["matnest"] = False
problem = fe.LinearVariationalProblem(a, L, w, bcs=[bc0, bc1])
solver = fe.LinearVariationalSolver(
problem, solver_parameters={'ksp_type':'preonly', 'pc_type':'lu'})
solver.solve()
This prints
Using firedrake-0.13.0+2351.gb7b1294 Traceback (most recent call last): File "mixed_poisson.py", line 46, in
solver.solve() File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 239, in solve work.copy(u) File "/usr/lib/python3.5/contextlib.py", line 77, in exit self.gen.throw(type, value, traceback) File "/home/zimmerman/firedrake/src/firedrake/firedrake/petsc.py", line 199, in inserted_options yield File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 239, in solve work.copy(u) File "/usr/lib/python3.5/contextlib.py", line 77, in exit self.gen.throw(type, value, traceback) File "/home/zimmerman/firedrake/src/PyOP2/pyop2/petsc_base.py", line 410, in vecscatter yield self._vec File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 238, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 538, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4396 in /tmp/pip-req-build-4g28hq0f/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 38 in /tmp/pip-req-build-4g28hq0f/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 723 in /tmp/pip-req-build-4g28hq0f/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-4g28hq0f/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-4g28hq0f/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 84 in /tmp/pip-req-build-4g28hq0f/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-4g28hq0f/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-4g28hq0f/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
So with the latest Firedrake, how does one disable matrix nesting, or generally how does one use a direct solver with mixed FE problems?
Dear Alex,
For mixed problems with LU you need to use the PETSc option "mat_type":"aij" to tell PETSc to make one big matrix instead of storing the blocks separately.
all the best
--cjc
From: Alexander G. Zimmerman notifications@github.com Sent: 13 November 2018 08:54:12 To: firedrakeproject/firedrake Cc: Subscribed Subject: Re: [firedrakeproject/firedrake] Direct solve fails for a mixed system (#431)
I'm trying to use a direct solver with a mixed problem and the solution to this issue no longer seems correct.
If I try to pass firedrake.LinearVariationalProblem.init does not have a nest keyword argument. Neither does firedrake.LinearVariationalSolver.solve. The same is true for the nonlinear counterparts. Setting fe.parameters["matnest"] = False does not change the behavior.
Here's a MWE based on the mixed Poisson tutorialhttps://www.firedrakeproject.org/demos/poisson_mixed.py.html:
import firedrake as fe
print("Using firedrake-" + fe.version)
mesh = fe.UnitSquareMesh(3, 3)
BDM = fe.FunctionSpace(mesh, "BDM", 1)
DG = fe.FunctionSpace(mesh, "DG", 0)
W = BDM * DG
sigma, u = fe.TrialFunctions(W)
tau, v = fe.TestFunctions(W)
x, y = fe.SpatialCoordinate(mesh)
f = fe.Function(DG).interpolate( 10*fe.exp(-(pow(x - 0.5, 2) + pow(y - 0.5, 2)) / 0.02))
dot, div, dx = fe.dot, fe.div, fe.dx
a = (dot(sigma, tau) + div(tau)u + div(sigma)v)*dx
L = - fvdx
bc0 = fe.DirichletBC(W.sub(0), fe.as_vector([0.0, -fe.sin(5*x)]), 1)
bc1 = fe.DirichletBC(W.sub(0), fe.as_vector([0.0, fe.sin(5*y)]), 2)
w = fe.Function(W)
fe.parameters["matnest"] = False
problem = fe.LinearVariationalProblem(a, L, w, bcs=[bc0, bc1])
solver = fe.LinearVariationalSolver( problem, solver_parameters={'ksp_type':'preonly', 'pc_type':'lu'})
solver.solve()
This prints
Using firedrake-0.13.0+2351.gb7b1294 Traceback (most recent call last): File "mixed_poisson.py", line 46, in solver.solve() File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 239, in solve work.copy(u) File "/usr/lib/python3.5/contextlib.py", line 77, in exit self.gen.throw(type, value, traceback) File "/home/zimmerman/firedrake/src/firedrake/firedrake/petsc.py", line 199, in inserted_options yield File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 239, in solve work.copy(u) File "/usr/lib/python3.5/contextlib.py", line 77, in exit self.gen.throw(type, value, traceback) File "/home/zimmerman/firedrake/src/PyOP2/pyop2/petsc_base.py", line 410, in vecscatter yield self._vec File "/home/zimmerman/firedrake/src/firedrake/firedrake/variational_solver.py", line 238, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 538, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4396 in /tmp/pip-req-build-4g28hq0f/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 38 in /tmp/pip-req-build-4g28hq0f/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 723 in /tmp/pip-req-build-4g28hq0f/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-4g28hq0f/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-4g28hq0f/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 84 in /tmp/pip-req-build-4g28hq0f/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-4g28hq0f/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-4g28hq0f/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
So with the latest Firedrake, how does one disable matrix nesting, or generally how does one use a direct solver with mixed FE problems?
— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/firedrakeproject/firedrake/issues/431#issuecomment-438186006, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AF0q8rChAHgt05qNo1BJRyct4xMW3kxZks5uuoi0gaJpZM4DTzsg.
@colinjcotter , thanks! There was only one more thing missing, which was setting "pc_factor_mat_solver_type": "mumps". So it works if I instantiate the solver with
solver = fe.LinearVariationalSolver(
problem,
solver_parameters = {
"ksp_type": "preonly",
"pc_type": "lu",
"mat_type": "aij",
"pc_factor_mat_solver_type": "mumps"})
Being able to search for "aij", I found this in the Navier-Stokes tutorial.
I'm trying to use a direct solver on a mixed system:
but I'm getting the following PETSc error:
Do I need to pass in some additional options when building PETSc, or is this unsupported in Firedrake? If it's the latter, please tag this as a 'wishlist' item.