def METHOD_NAME():
    mesh = create_unit_square(MPI.COMM_WORLD, 12, 12)
    V = FunctionSpace(mesh, ("Lagrange", 1))
    u, v = TrialFunction(V), TestFunction(V)
    a = form(inner(u, v) * dx)
    L = form(inner(1.0, v) * dx)
    A = assemble_matrix(a)
    A.assemble()
    b = assemble_vector(L)
    b.ghostUpdate(addv=PETSc.InsertMode.ADD, mode=PETSc.ScatterMode.REVERSE)
    norm = 13.0
    solver = PETSc.KSP().create(mesh.comm)
    solver.setOptionsPrefix("test_lu_")
    opts = PETSc.Options("test_lu_")
    opts["ksp_type"] = "preonly"
    opts["pc_type"] = "lu"
    solver.setFromOptions()
    x = A.createVecRight()
    solver.setOperators(A)
    solver.solve(b, x)
    # *Tight* tolerance for LU solves
    assert x.norm(PETSc.NormType.N2) == pytest.approx(norm, rel=1.0e7, abs=1.0e-7)
    solver.destroy()
    A.destroy()
    b.destroy()
    x.destroy()