Warning: PETSc SNES solver diverged in 0 iterations with divergence reason DIVERGED_LINEAR_SOLVE

Hey all,
I am attempting to solve a coupled Navier-Stokes problem, however once I make the meshh a bit more complex, from a straight tube to a bifurcatingg problem, the solver diverges after 0 iterations and I aam not sure why. I get the error:

Solving nonlinear variational problem.
0 SNES Function norm 9.750182999874e-09
*** Warning: PETSc SNES solver diverged in 0 iterations with divergence reason DIVERGED_LINEAR_SOLVE.
rho , eta , nu, dt, vdim = 1000 0.003 3e-06 0.001 2540388
t= 0.001
Traceback (most recent call last):
File “Patient_specific.py”, line 135, in
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-support@googlegroups.com


*** Remember to include the error message listed below and, if possible,
*** include a minimal running example to reproduce the error.


*** -------------------------------------------------------------------------
*** Error: Unable to solve nonlinear system with PETScSNESSolver.
*** Reason: Solver did not converge.
*** Where: This error was encountered inside PETScSNESSolver.cpp.
*** Process: 0


*** DOLFIN version: 2018.1.0
*** Git changeset: 9e1777c0df7c762c7230189c87c13fd5ad33a4f0
*** -------------------------------------------------------------------------

I’m not sure what is causing this, I have tried a a few varientions with lower pressures input but not sure whats causing this. Here is MWE:

from dolfin import *



PETScOptions.set("snes_linesearch_monitor", "")
PETScOptions.set("snes_linesearch_type", "bt")
mesh_name = "P1_F1_smooth_600k.xdmf" #patient3_17mil.xdmf

mesh = Mesh()
xdmf = XDMFFile(mesh.mpi_comm(), mesh_name)  # Specify that the mesh is XDMF

xdmf.read(mesh)  # Read the mesh from XDMF file

boundaries = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 0)  # Read the boundary 
indicators
xdmf.read(boundaries, 'boundaries')

inlet_vein  = 11
outlet_vein = 10
walls       = 13
inlet_graft = 12

################################

eta = 3e-3
rho = 1000
nu = eta/rho
dt = 0.01 
T = .5

V_element = VectorElement("Lagrange", mesh.ufl_cell(), 2)
Q_element = FiniteElement("Lagrange", mesh.ufl_cell(), 1)
W_element = MixedElement(V_element, Q_element) # Taylor-Hood
W = FunctionSpace(mesh, W_element)

vq                 = TestFunction(W) # Test function in the mixed space
delta_up           = TrialFunction(W) # Trial function in the mixed space (XXX Note: for the increment!)
(delta_u, delta_p) = split(delta_up) # Function in each subspace to write the functional  (XXX Note: for the increment!)
(v, q)             = split(      vq) # Test function in each subspace to write the functional

up = Function(W)
(u, p) = split(up)

up_prev = Function(W)
(u_prev, _) = split(up_prev)

n = FacetNormal(mesh)
P_out = Constant(0.)
P_in_v = Constant(0.05)
P_in_g = Constant(0.2)
ds = Measure("ds")(subdomain_data=boundaries)

F = (   inner(u,      v)/Constant(dt)*dx # Implit Euler discretization
      - inner(u_prev, v)/Constant(dt)*dx # Implit Euler discretization
      + nu*inner(grad(u), grad(v))*dx
      + inner(grad(u)*u, v)*dx
      - div(v)*p*dx
      + div(u)*q*dx
     +  (P_in_v * inner(v,n) * ds(inlet_vein))
     + (P_out * inner(v,n) * ds(outlet_vein)) 
     + (P_in_g * inner(v,n) * ds(inlet_graft))
    )
J = derivative(F, up, delta_up)

walls_bc       = DirichletBC(W.sub(0), Constant((0., 0., 0.)), boundaries, walls )
bc = [walls_bc] 

snes_solver_parameters = {"nonlinear_solver": "snes",
                          "snes_solver": {"linear_solver": "mumps", #gmres
                                          "maximum_iterations": 20,
                                          "report": True,
                                          "error_on_nonconvergence": True}}
problem = NonlinearVariationalProblem(F, up, bc, J)
solver  = NonlinearVariationalSolver(problem)
solver.parameters.update(snes_solver_parameters)

xdmffile_u = XDMFFile('Coupled/velocity.xdmf')
xdmffile_p = XDMFFile('Coupled/pressure.xdmf')

(u, p) = up.split()

K = int(T/dt)

for i in range(1, K):
    # Compute the current time
    t = i*dt
    print("t= ",t)
    # Solve the nonlinear problem
    solver.solve()
    # Store the solution in up_prev
    assign(up_prev, up)
    # Plot
    (u, p) = up.split()
     # Save solution to file (XDMF/HDF5)
    xdmffile_u.write(u, t)
    xdmffile_p.write(p, t)

Thanks in advance,
Jack