@@ -0,0 +1,1803 @@
+from firedrake import *%0Afrom firedrake.slate.preconditioners import create_schur_nullspace%0Aimport numpy as np%0A%0Amesh = UnitCubedSphereMesh(2)%0Amesh.init_cell_orientations(SpatialCoordinate(mesh))%0A%0An = FacetNormal(mesh)%0A%0AV = FunctionSpace(mesh, %22RTCF%22, 1)%0AQ = FunctionSpace(mesh, %22DG%22, 0)%0A%0AW = V*Q%0Asigma, u = TrialFunctions(W)%0Atau, v = TestFunctions(W)%0Aa = (inner(sigma, tau) + div(sigma)*v + div(tau)*u)*dx%0A%0AW_d = FunctionSpace(mesh,%0A                    MixedElement(%5BBrokenElement(Vi.ufl_element())%0A                                  for Vi in W%5D))%0A%0Aatilde = Tensor(replace(a, dict(zip(a.arguments(),%0A                                    (TestFunction(W_d),%0A                                     TrialFunction(W_d))))))%0A%0AVt = FunctionSpace(mesh, %22HDiv Trace%22, 0)%0Agamma = TestFunction(Vt)%0A%0Asigma, _ = TrialFunctions(W_d)%0A%0AK = Tensor(gamma('+') * dot(sigma, n) * dS)%0A%0AA = assemble(a, mat_type=%22aij%22)%0Anullspace = MixedVectorSpaceBasis(W, %5BW%5B0%5D, VectorSpaceBasis(constant=True)%5D)%0Anullspace._build_monolithic_basis()%0AA.petscmat.setNullSpace(nullspace._nullspace)%0A%0ASnullsp = create_schur_nullspace(A.petscmat, K * atilde.inv,%0A                                 W, W_d, Vt,%0A                                 COMM_WORLD)%0A%0Av = Snullsp.getVecs()%5B0%5D%0Aprint %22Computed nullspace of S (min, max, norm)%22, v.array_r.min(), v.array_r.max(), v.norm()%0A%0AS = K * atilde.inv * K.T%0A%0Au, s, v = np.linalg.svd(assemble(S, mat_type=%22aij%22).M.values)%0A%0Asingular_vector = v%5B-1%5D%0Aprint %22Actual nullspace of S (min, max, norm)%22, singular_vector.min(), singular_vector.max(), np.linalg.norm(singular_vector)%0A%0Au, s, v = np.linalg.svd(A.M.handle%5B:, :%5D)%0A%0Aoffset = V.dof_dset.size%0A%0Asingular_vector = v%5B-1%5D%5Boffset:%5D%0Aprint %22Nullspace of original operator (min, max, norm)%22, singular_vector.min(), singular_vector.max(), np.linalg.norm(singular_vector)%0A

Add test for computing nullspace of the Schur operator