firedrakeproject / firedrake

Firedrake is an automated system for the portable solution of partial differential equations using the finite element method (FEM)
https://firedrakeproject.org
Other
482 stars 157 forks source link

BUG: Arity zero `Action` weights ignored in assembly #3597

Closed jrmaddison closed 4 weeks ago

jrmaddison commented 1 month ago

Describe the bug Assembling a multiple of an arity zero Action discards the weight, and returns the unscaled action.

Steps to Reproduce

from firedrake import *

mesh = UnitIntervalMesh(1)
space = FunctionSpace(mesh, "Discontinuous Lagrange", 0)
test = TestFunction(space)
trial = TrialFunction(space)

u = Function(space, name="u").interpolate(Constant(1.0))
v = Cofunction(space.dual(), name="v")
assemble(test * dx, tensor=v)

for f in (1.0, 0.5, -1.0):
    print(f"{(f, assemble(f * v(u)))=}")

displays

(f, assemble(f * v(u)))=(1.0, 1.0)
(f, assemble(f * v(u)))=(0.5, 1.0)
(f, assemble(f * v(u)))=(-1.0, 1.0)

Expected behavior The result should be scaled by the weight.

Error message No error.

Environment: Ubuntu 22.04, recent Firedrake build.

jrmaddison commented 4 weeks ago

I think the weight is discarded here

https://github.com/firedrakeproject/firedrake/blob/b6af9b2f52b8d37118b6f7e57efebc1500c5a432/firedrake/assemble.py#L636-L638

For a FormSum expr.ufl_operands excludes the weight (expr.weights()).

jrmaddison commented 4 weeks ago

Sorry, that's wrong, it's from here.

https://github.com/firedrakeproject/firedrake/blob/b6af9b2f52b8d37118b6f7e57efebc1500c5a432/firedrake/assemble.py#L466-L467

The sum should include the expr.weights().