I am trying to solve the navier stokes equations using an implicit method combined with newtons method to deal with the nonlinearity. When I try to do so I get the error seen at the bottom of this question. I've included a simplified version of my code below which recreates the error.
from dolfin import *
from mshr import *
import numpy as np
circle_outx = 0.0
circle_outy = 0.0
circle_outr = 1.0
circle_inx = 0.5
circle_iny = 0.0
circle_inr = 0.1
domain = Circle(Point(circle_outx,circle_outy),circle_outr) - Circle(Point(circle_inx,circle_iny),circle_inr)
mesh = generate_mesh ( domain, 30 )
#Define Function Spaces
V = VectorFunctionSpace(mesh, "CG", 2)
Q = FunctionSpace(mesh, "CG", 1)
W = V * Q
# Boundaries
def u0_boundary(x, on_boundary):
return on_boundary
noslip = Constant((0.0, 0.0))
bc0 = DirichletBC(W.sub(0),noslip,u0_boundary)
bc1 = DirichletBC(W.sub(1),0.0,u0_boundary)
bcs = [bc0, bc1]
#Functions and TestFunctions
u = Function(V)
p = Function(Q)
u0 = Function(V)
p0 = Function(Q)
w0 = Function(W)
w = Function(W)
(v, q) = TestFunctions(W)
# functions for stokes
(us, ps) = TrialFunctions(W)
(vs, qs) = TestFunctions(W)
# weak form Stokes equation
fstokes = Expression(("-1*4*x[1]*(1 - x[0]*x[0] - x[1]*x[1])","4*x[0]*(1 - x[0]*x[0] - x[1]*x[1])"))
Stokes = (inner(grad(us), grad(vs)) - div(vs)*ps + qs*div(us))*dx
LStokes = inner(fstokes, vs)*dx
#Get the initial condition
solve(Stokes == LStokes, w0, bcs,solver_parameters=dict(linear_solver="lu"))
(uinit, pinit) = w0.split(True)
u0.assign(uinit)
p0.assign(pinit)
#This is the nonlinear equation I can't solve
F = inner(u,v)*dx + inner(grad(u)*u,v)*dx - inner(u0,v) *dx
J = derivative ( F, w )
solve ( F == 0, w, bcs, J = J )
Solving nonlinear variational problem.
Newton iteration 0: r (abs) = 2.292e-03 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-09)
[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Matrix is missing diagonal entry in row 1!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.3.0, Patch 3, Wed Aug 29 11:26:24 CDT 2012
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Unknown Name on a darwin14. named Michaels-MacBook-Pro-3.local by michaelschneier Fri Oct 16 17:09:28 2015
[0]PETSC ERROR: Libraries linked from /Users/johannr/fenics-1.5.0/local/lib
[0]PETSC ERROR: Configure run at Mon Jan 12 22:17:07 2015
[0]PETSC ERROR: Configure options --prefix=/Users/johannr/fenics-1.5.0/local COPTFLAGS=-O2 --with-debugging=0 --with-clanguage=cxx --with-c-support=1 --with-blas-lapack-dir=/usr --with-umfpack=1 --with-umfpack-include=/Users/johannr/fenics-1.5.0/local/include/suitesparse --with-umfpack-lib="[/Users/johannr/fenics-1.5.0/local/lib/libumfpack.a,/Users/johannr/fenics-1.5.0/local/lib/libamd.a]" --with-spooles=1 --with-spooles-include=/Users/johannr/fenics-1.5.0/local/include --with-spooles-lib=/Users/johannr/fenics-1.5.0/local/lib/libspooles.a --with-ptscotch=1 --with-ptscotch-dir=/Users/johannr/fenics-1.5.0/local --with-ml=1 --with-ml-include=/Users/johannr/fenics-1.5.0/local/include/trilinos --with-ml-lib=/Users/johannr/fenics-1.5.0/local/lib/libml.dylib --with-hdf5=1 --with-hdf5-dir=/Users/johannr/fenics-1.5.0/local --with-x=0 -with-x11=0 --with-fortran=0 --with-shared-libraries=1 PETSC_DIR=/Users/johannr/fenics-1.5.0/fenics-superbuild/build-fenics/CMakeExternals/src/PETSc PETSC_ARCH=darwin14.0.0-cxx-opt
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: MatZeroRows_SeqAIJ() line 1556 in /Users/johannr/fenics-1.5.0/fenics-superbuild/build-fenics/CMakeExternals/src/PETSc/src/mat/impls/aij/seq/aij.c
[0]PETSC ERROR: MatZeroRowsLocal() line 5652 in /Users/johannr/fenics-1.5.0/fenics-superbuild/build-fenics/CMakeExternals/src/PETSc/src/mat/interface/matrix.c
Traceback (most recent call last):
File "testing.py", line 56, in
solve ( F == 0, w, bcs, J = J )
File "/Applications/FEniCS.app/Contents/Resources/lib/python2.7/site-packages/dolfin/fem/solving.py", line 269, in solve
_solve_varproblem(*args, **kwargs)
File "/Applications/FEniCS.app/Contents/Resources/lib/python2.7/site-packages/dolfin/fem/solving.py", line 317, in _solve_varproblem
solver.solve()
RuntimeError:
*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
*** fenics@fenicsproject.org
*** Remember to include the error message listed below and, if possible,
*** include a minimal running example to reproduce the error.
*** -------------------------------------------------------------------------
*** Error: Unable to set given (local) rows to identity matrix.
*** Reason: some diagonal elements not preallocated (try assembler option keep_diagonal).
*** Where: This error was encountered inside PETScMatrix.cpp.
*** Process: unknown
*** DOLFIN version: 1.5.0
*** Git changeset: f467b66dcfd821ec20e9f9070c7cef5a991dbc42
*** -----------------------------------------------------
