BUG: Solving the Stokes equations defined with TrialFunction returns an error

[diego-hayashi]

Describe the bug When solving the Stokes equations that were defined with TrialFunction, Firedrake returns an error.

Steps to Reproduce Run the following code, which is based on stokes-topology-rol-firedrake.py, and set to run with TrialFunction:

from firedrake import *

N = 20
mesh = RectangleMesh(N, N, 1.5, 1, diagonal = "right")

U_h = VectorElement("CG", mesh.ufl_cell(), 2)
P_h = FiniteElement("CG", mesh.ufl_cell(), 1)
W = FunctionSpace(mesh, U_h*P_h)

(x, y) = SpatialCoordinate(mesh)
l = 1.0/6.0
gbar = 1.0
cond1 = And(gt(y, (1.0/4 - l/2)), lt(y, (1.0/4 + l/2)))
val1  = gbar*(1 - (2*(y-0.25)/l)**2)
cond2 = And(gt(y, (3.0/4 - l/2)), lt(y, (3.0/4 + l/2)))
val2  = gbar*(1 - (2*(y-0.75)/l)**2)
inflow_outflow = conditional(cond1, val1, conditional(cond2, val2, 0))

w = Function(W)
(u, p) = TrialFunctions(W)
(v, q) = TestFunctions(W)

F = inner(grad(u), grad(v)) * dx + inner(grad(p), v) * dx  + inner(div(u), q) * dx
bcs = [DirichletBC(W.sub(0).sub(1), 0, "on_boundary"), DirichletBC(W.sub(0).sub(0), inflow_outflow, "on_boundary")]

a = lhs(F)
L = rhs(F)
solve(a == L, w, bcs = bcs)

Expected behavior The code should have run fine.

Error message

Traceback (most recent call last):
  File "test.py", line 36, in <module>
    solve(a == L, w, bcs = bcs)
    ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "petsc4py/PETSc/Log.pyx", line 188, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "petsc4py/PETSc/Log.pyx", line 189, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "firedrake/src/firedrake/firedrake/adjoint_utils/solving.py", line 57, in wrapper
    output = solve(*args, **kwargs)
             ^^^^^^^^^^^^^^^^^^^^^^
  File "firedrake/src/firedrake/firedrake/solving.py", line 140, in solve
    _solve_varproblem(*args, **kwargs)
  File "firedrake/src/firedrake/firedrake/solving.py", line 167, in _solve_varproblem
    problem = vs.LinearVariationalProblem(eq.lhs, eq.rhs, u, bcs, Jp,
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "petsc4py/PETSc/Log.pyx", line 188, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "petsc4py/PETSc/Log.pyx", line 189, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "firedrake/src/firedrake/firedrake/variational_solver.py", line 370, in __init__
    raise ValueError("Provided RHS is not a linear form")

Additional Info If the line L = rhs(F), which is a form of the empty type, is changed to L = Constant(0.)*q*dx, Firedrake is able to solve the equation.

[diego-hayashi]

One possible fix:

diff --git a/firedrake/variational_solver.py b/firedrake/variational_solver.py
index bba281bc1..430a9c311 100644
--- a/firedrake/variational_solver.py
+++ b/firedrake/variational_solver.py
@@ -366,7 +366,7 @@ class LinearVariationalProblem(NonlinearVariationalProblem):
         else:
             if not isinstance(L, (ufl.BaseForm, slate.slate.TensorBase)):
                 raise TypeError("Provided RHS is a '%s', not a Form or Slate Tensor" % type(L).__name__)
-            if len(L.arguments()) != 1:
+            if len(L.arguments()) > 1:
                 raise ValueError("Provided RHS is not a linear form")
             F = ufl_expr.action(J, u) - L