In case anyone has the same problem, I had a similar issue and solved it on my end.
On the subject of boundary conditions, my understanding of the problem is that under the hood, when a DirichletBC is specified at assembly, dolfin cancels the rows and columns associated with the relevant degrees of freedom, then set its diagonal terms to 1. In other words, given a boundary \Gamma it modifies A so that Ax=b also imposes x_\Gamma=b_\Gamma
Another very important physical thing about the code above is the lack of a mass matrix. Physically, you have the Navier-Stokes equations as :
\begin{cases}\partial_tu+(u\cdot \nabla)u=-\nabla p+\frac{1}{Re}\Delta u\\ \nabla\cdot u=0\end{cases}\Leftrightarrow\left[\begin{array}{cc}I&0\\0&0\end{array}\right]\partial_t\left[\begin{array}{c}u\\p\end{array}\right]=\mathcal{F}(w)\Leftrightarrow\mathcal{N}\partial_tw=\mathcal{F}(w)
Indeed one can infer from your weak form that your equations are incompressible. This means your eigenvalue problem really writes as a general eigenvalue problem :
\sigma\mathcal{N}\hat{w}=\mathcal{J}(\bar{w})\hat{w}
This \mathcal{N} matrix is fairly trivial to write in fenics, with the specificity that if you want your boundary conditions to be specified. You’re studying perturbations so your boundary conditions are homogeneous and you want \mathcal{N} to cancel on your boundary \Gamma.
I’m coding in dolfinx now, but I think in dolfin the most elegant way to fix your code would be to replace the central part by :
A = PETScMatrix()
B = PETScMatrix()
#Assembling the jacobian of the weak form evaluated at the baseflow previously calculated and stored in w
J = derivative(F, w)
really_meaningful = inner(u,v)*dx
assemble_system(J, really_meaningful, bcs, A_tensor = A, B_tensor=B)
#Solving the standard eigenvalue problem, looking for the largest real part
eigensolver = SLEPcEigenSolver(-A,B)
Hope that is of some use to someone someday.