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
483 stars 157 forks source link

Riesz representative changed from l2 to L2 #3576

Closed APaganini closed 1 month ago

APaganini commented 1 month ago

Description

Ensure that default Riesz representative is L2 and not l2 by changing line 228 of firedrake/adjoint_utils/function.py from

riesz_representation = options.get("riesz_representation", "l2") to

riesz_representation = options.get("riesz_representation", "L2") This fixes #3575

Ig-dolci commented 1 month ago

Thank you for catching this bug! Perhaps a test can be added to this PR? I just wrote this:

    # Make sure the derivative provides the right computation 
    # for the Riesz representation with the H1, L2, l2 and None options.
    mesh = UnitIntervalMesh(1)
    space = FunctionSpace(mesh, "CG", 1)
    f = Function(space)
    x = SpatialCoordinate(mesh)
    f.interpolate(x[0])
    J = assemble((f ** 2) * dx)
    with stop_annotating():
        v = TestFunction(space)
        u = TrialFunction(space)
        dJdu_cofunction = assemble(2 * inner(f, v) * dx)

        # Riesz representation with l2
        dJdu_function_l2 = Function(space, val=dJdu_cofunction.dat)

        # Riesz representation with H1
        a = firedrake.inner(u, v)*firedrake.dx \
            + firedrake.inner(firedrake.grad(u), firedrake.grad(v))*firedrake.dx
        dJdu_function_H1 = Function(space)
        solve(a == dJdu_cofunction, dJdu_function_H1)

        # Riesz representation with L2
        a = firedrake.inner(u, v)*firedrake.dx
        dJdu_function_L2 = Function(space)
        solve(a == dJdu_cofunction, dJdu_function_L2)

    rf = ReducedFunctional(J, Control(f))
    assert np.allclose(
        rf.derivative(options={"riesz_representation": None}).dat.data_ro, dJdu_cofunction.dat.data_ro
        )
    assert np.allclose(
        rf.derivative(options={"riesz_representation": "H1"}).dat.data_ro, dJdu_function_H1.dat.data_ro
        )
    assert np.allclose(
        rf.derivative(options={"riesz_representation": "l2"}).dat.data_ro, dJdu_function_l2.dat.data_ro
        )
    # L2 is the default option
    assert np.allclose(rf.derivative().dat.data_ro, dJdu_function_L2.dat.data_ro)
Ig-dolci commented 1 month ago

This PR will be closed as PR 3579 addresses the same changes.