From 4e404cc6f3ed0b1761e3e16be7caad1736b4b0e4 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Sun, 18 May 2025 17:19:39 +0100
Subject: [PATCH 01/11] Re-enable adjoint interpolation

---
 animate/interpolation.py   | 4 +---
 test/test_interpolation.py | 8 ++++----
 2 files changed, 5 insertions(+), 7 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index d2e1102a..7d233567 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -254,9 +254,7 @@ def _transfer_adjoint(target_b, source_b, transfer_method, **kwargs):
     # Apply adjoint transfer operator to each component
     for i, (t_b, s_b) in enumerate(zip(target_b_split, source_b_split)):
         if transfer_method == "interpolate":
-            raise NotImplementedError(
-                "Adjoint of interpolation operator not implemented."
-            )  # TODO (#113)
+            s_b.interpolate(t_b, adjoint=True, **kwargs)
         elif transfer_method == "project":
             ksp = petsc4py.KSP().create()
             ksp.setOperators(assemble_mass_matrix(t_b.function_space(), lumped=bounded))
diff --git a/test/test_interpolation.py b/test/test_interpolation.py
index 1cdfd65a..485937c8 100644
--- a/test/test_interpolation.py
+++ b/test/test_interpolation.py
@@ -255,7 +255,7 @@ def test_transfer_same_space(self, transfer_method):
         expected = source
         self.assertAlmostEqual(errornorm(expected, target), 0)
 
-    @parameterized.expand(["project"])  # TODO: interpolate (#113)
+    @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_space_adjoint(self, transfer_method):
         pytest.skip()  # TODO: (#114)
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
@@ -279,7 +279,7 @@ def test_transfer_same_space_mixed(self, transfer_method):
         expected = source
         self.assertAlmostEqual(errornorm(expected, target), 0)
 
-    @parameterized.expand(["project"])  # TODO: interpolate (#113)
+    @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_space_mixed_adjoint(self, transfer_method):
         pytest.skip()  # TODO: (#114)
         P1 = FunctionSpace(self.source_mesh, "CG", 1)
@@ -307,7 +307,7 @@ def test_transfer_same_mesh(self, transfer_method):
             expected = Function(Vt).project(source)
         self.assertAlmostEqual(errornorm(expected, target), 0)
 
-    @parameterized.expand(["project"])  # TODO: interpolate (#113)
+    @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_mesh_adjoint(self, transfer_method):
         pytest.skip()  # TODO: (#114)
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
@@ -343,7 +343,7 @@ def test_transfer_same_mesh_mixed(self, transfer_method):
             e2.project(s2)
         self.assertAlmostEqual(errornorm(expected, target), 0)
 
-    @parameterized.expand(["project"])  # TODO: interpolate (#113)
+    @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_mesh_mixed_adjoint(self, transfer_method):
         pytest.skip()  # TODO: (#114)
         P1 = FunctionSpace(self.source_mesh, "CG", 1)

From 78b9efea836ca443e789226d402f8facb683700d Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Sun, 18 May 2025 17:49:22 +0100
Subject: [PATCH 02/11] Accept rvalue for (co)function2(co)function

---
 animate/utility.py | 14 ++++++++++----
 1 file changed, 10 insertions(+), 4 deletions(-)

diff --git a/animate/utility.py b/animate/utility.py
index 16d450a6..1ac0c629 100644
--- a/animate/utility.py
+++ b/animate/utility.py
@@ -279,14 +279,17 @@ def assemble_mass_matrix(space, norm_type="L2", lumped=False):
         return mass_matrix.createDiagonal(x)
 
 
-def cofunction2function(cofunc):
+def cofunction2function(cofunc, func=None):
     """
     :arg cofunc: a cofunction
     :type cofunc: :class:`firedrake.cofunction.Cofunction`
+    :kwarg func: a function for the return value
+    :type func: :class:`firedrake.function.Function`
     :returns: a function with the same underyling data
     :rtype: :class:`firedrake.function.Function`
     """
-    func = ffunc.Function(cofunc.function_space().dual())
+    if func is None:
+        func = ffunc.Function(cofunc.function_space().dual())
     if isinstance(func.dat.data_with_halos, tuple):
         for i, arr in enumerate(func.dat.data_with_halos):
             arr[:] = cofunc.dat.data_with_halos[i]
@@ -295,14 +298,17 @@ def cofunction2function(cofunc):
     return func
 
 
-def function2cofunction(func):
+def function2cofunction(func, cofunc=None):
     """
     :arg func: a function
     :type func: :class:`firedrake.function.Function`
+    :kwarg cofunc: a cofunction for the return value
+    :type cofunc: :class:`firedrake.cofunction.Cofunction`
     :returns: a cofunction with the same underlying data
     :rtype: :class:`firedrake.cofunction.Cofunction`
     """
-    cofunc = firedrake.Cofunction(func.function_space().dual())
+    if cofunc is None:
+        cofunc = firedrake.Cofunction(func.function_space().dual())
     if isinstance(cofunc.dat.data_with_halos, tuple):
         for i, arr in enumerate(cofunc.dat.data_with_halos):
             arr[:] = func.dat.data_with_halos[i]

From 9148fe4baa61fb0c25bc37d1cbf3bda31dee270f Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:28:59 +0100
Subject: [PATCH 03/11] Rework interpolation

---
 animate/interpolation.py | 243 ++++++++++++---------------------------
 1 file changed, 75 insertions(+), 168 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index 7d233567..7b4475ce 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -66,63 +66,55 @@ def transfer(source, target_space, transfer_method="project", **kwargs):
         raise TypeError(
             "Second argument must be a FunctionSpace, Function, or Cofunction."
         )
-    if isinstance(source, firedrake.Cofunction):
-        return _transfer_adjoint(source, target, transfer_method, **kwargs)
-    elif source.function_space() == target.function_space():
-        return target.assign(source)
+    if transfer_method == "interpolate":
+        return interpolate(source, target, **kwargs)
     else:
-        return _transfer_forward(source, target, transfer_method, **kwargs)
+        return project(source, target, **kwargs)
 
 
-@PETSc.Log.EventDecorator()
-def interpolate(source, target_space, **kwargs):
-    r"""
-    Overload function :func:`firedrake.__future__.interpolate` to account for the case
-    of two mixed function spaces defined on different meshes and for the adjoint
-    interpolation operator when applied to :class:`firedrake.cofunction.Cofunction`\s.
-
-    :arg source: the function to be transferred
-    :type source: :class:`firedrake.function.Function` or
-        :class:`firedrake.cofunction.Cofunction`
-    :arg target_space: the function space which we seek to transfer onto, or the
-        function or cofunction to use as the target
-    :type target_space: :class:`firedrake.functionspaceimpl.FunctionSpace`,
-        :class:`firedrake.function.Function` or :class:`firedrake.cofunction.Cofunction`
-    :returns: the transferred function
-    :rtype: :class:`firedrake.function.Function` or
-        :class:`firedrake.cofunction.Cofunction`
-
-    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate`
-    """
-    return transfer(source, target_space, transfer_method="interpolate", **kwargs)
+def _validate_consistent_spaces(Vs, Vt):
+    if hasattr(Vs, "num_sub_spaces"):
+        if not hasattr(Vt, "num_sub_spaces"):
+            raise ValueError(
+                "Source space has multiple components but target space does not."
+            )
+        if Vs.num_sub_spaces() != Vt.num_sub_spaces():
+            raise ValueError(
+                "Inconsistent numbers of components in source and target spaces:"
+                f" {Vs.num_sub_spaces()} vs. {Vt.num_sub_spaces()}."
+            )
+    elif hasattr(Vt, "num_sub_spaces"):
+        raise ValueError(
+            "Target space has multiple components but source space does not."
+        )
+    if Vs._dual and not Vt._dual:
+        raise ValueError("Source space is a dual space but target space is not.")
+    if not Vs._dual and Vt._dual:
+        raise ValueError("Target space is a dual space but source space is not.")
 
 
 @PETSc.Log.EventDecorator()
-def project(source, target_space, **kwargs):
+def interpolate(source, target, **kwargs):
     r"""
-    Overload function :func:`firedrake.projection.project` to account for the case of
-    two mixed function spaces defined on different meshes and for the adjoint
-    projection operator when applied to :class:`firedrake.cofunction.Cofunction`\s.
+    Overload :func:`firedrake.__future__.interpolate` to account for the case of mixed
+    function spaces.
 
-    For details on the approach for achieving boundedness through mass lumping and
-    post-processing, see :cite:`Farrell:2009`.
-
-    :arg source: the function to be transferred
+    :arg source: the function or cofunction to be transferred
     :type source: :class:`firedrake.function.Function` or
         :class:`firedrake.cofunction.Cofunction`
-    :arg target_space: the function space which we seek to transfer onto, or the
-        function or cofunction to use as the target
-    :type target_space: :class:`firedrake.functionspaceimpl.FunctionSpace`,
-        :class:`firedrake.function.Function` or :class:`firedrake.cofunction.Cofunction`
-    :returns: the transferred function
-    :rtype: :class:`firedrake.function.Function` or
+    :arg target: the function or cofunction to use as the target, which is modified in
+        place
+    :type target: :class:`firedrake.function.Function` or
         :class:`firedrake.cofunction.Cofunction`
-    :kwarg bounded: apply mass lumping to the mass matrix to ensure boundedness
-    :type bounded: :class:`bool`
 
-    Extra keyword arguments are passed to :func:`firedrake.projection.project`.
+    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate`
     """
-    return transfer(source, target_space, transfer_method="project", **kwargs)
+    _validate_consistent_spaces(source.function_space(), target.function_space())
+    if hasattr(target.function_space(), "num_sub_spaces"):
+        for s, t in zip(source.subfunctions, target.subfunctions):
+            t.interpolate(s, **kwargs)
+    else:
+        target.interpolate(source, **kwargs)
 
 
 # TODO: Reimplement by introducing a LumpedSupermeshProjector subclass of
@@ -148,146 +140,61 @@ def _supermesh_project(source, target, bounded=False):
 
 
 @PETSc.Log.EventDecorator()
-def _transfer_forward(source, target, transfer_method, **kwargs):
-    """
-    Apply mesh-to-mesh transfer operator to a Function.
+def project(source, target, **kwargs):
+    r"""
+    Overload :func:`firedrake.projection.project` to account for the case of mixed
+    function spaces.
 
-    This function extends the functionality of :func:`firedrake.__future__.interpolate`
-    and :func:`firedrake.projection.project` to account for mixed spaces.
+    For details on the approach for achieving boundedness through mass lumping and
+    post-processing, see :cite:`Farrell:2009`.
 
-    :arg source: the Function to be transferred
-    :type source: :class:`firedrake.function.Function`
-    :arg target: the Function which we seek to transfer onto
-    :type target: :class:`firedrake.function.Function`
-    :kwarg transfer_method: the method to use for the transfer. Options are
-        "interpolate" (default) and "project"
-    :type transfer_method: str
+    :arg source: the function or cofunction to be transferred
+    :type source: :class:`firedrake.function.Function` or
+        :class:`firedrake.cofunction.Cofunction`
+    :arg target: the function or cofunction to transfer onto, which is modified in
+        place
+    :type target: :class:`firedrake.function.Function` or
+        :class:`firedrake.cofunction.Cofunction`
     :kwarg bounded: apply mass lumping to the mass matrix to ensure boundedness
-        (project only)
     :type bounded: :class:`bool`
-    :returns: the transferred Function
-    :rtype: :class:`firedrake.function.Function`
 
-    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate` or
-        :func:`firedrake.projection.project`.
+    Extra keyword arguments are passed to :func:`firedrake.projection.project`.
     """
-    is_project = transfer_method == "project"
-    bounded = is_project and kwargs.pop("bounded", False)
-    Vs = source.function_space()
-    Vt = target.function_space()
-    _validate_matching_spaces(Vs, Vt)
-    assert isinstance(target, firedrake.Function)
-    if hasattr(Vt, "num_sub_spaces"):
-        for s, t in zip(source.subfunctions, target.subfunctions):
-            if transfer_method == "interpolate":
-                t.interpolate(s, **kwargs)
-            elif transfer_method == "project":
+    _validate_consistent_spaces(source.function_space(), target.function_space())
+    if not source.function_space()._dual:
+        if hasattr(target.function_space(), "num_sub_spaces"):
+            for s, t in zip(source.subfunctions, target.subfunctions):
                 if bounded:
                     _supermesh_project(s, t, bounded=True)
                 else:
                     t.project(s, **kwargs)
-            else:
-                raise ValueError(
-                    f"Invalid transfer method: {transfer_method}."
-                    " Options are 'interpolate' and 'project'."
-                )
-    else:
-        if transfer_method == "interpolate":
-            target.interpolate(source, **kwargs)
-        elif transfer_method == "project":
+        else:
             if bounded:
                 _supermesh_project(source, target, bounded=True)
             else:
                 target.project(source, **kwargs)
-        else:
-            raise ValueError(
-                f"Invalid transfer method: {transfer_method}."
-                " Options are 'interpolate' and 'project'."
-            )
-    return target
-
-
-@PETSc.Log.EventDecorator()
-def _transfer_adjoint(target_b, source_b, transfer_method, **kwargs):
-    """
-    Apply an adjoint mesh-to-mesh transfer operator to a Cofunction.
-
-    :arg target_b: seed Cofunction from the target space of the forward projection
-    :type target_b: :class:`firedrake.cofunction.Cofunction`
-    :arg source_b: output Cofunction from the source space of the forward projection
-    :type source_b: :class:`firedrake.cofunction.Cofunction`
-    :kwarg transfer_method: the method to use for the transfer. Options are
-        "interpolate" (default) and "project"
-    :type transfer_method: str
-    :kwarg bounded: apply mass lumping to the mass matrix to ensure boundedness
-        (project only)
-    :type bounded: :class:`bool`
-    :returns: the transferred Cofunction
-    :rtype: :class:`firedrake.cofunction.Cofunction`
-
-    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate` or
-        :func:`firedrake.projection.project`.
-    """
-    is_project = transfer_method == "project"
-    bounded = is_project and kwargs.pop("bounded", False)
-
-    # Map to Functions to apply the adjoint transfer
-    if not isinstance(target_b, firedrake.Function):
-        target_b = cofunction2function(target_b)
-    if not isinstance(source_b, firedrake.Function):
-        source_b = cofunction2function(source_b)
-
-    Vt = target_b.function_space()
-    Vs = source_b.function_space()
-    if Vs == Vt:
-        source_b.assign(target_b)
-        return function2cofunction(source_b)
-
-    _validate_matching_spaces(Vs, Vt)
-    if hasattr(Vs, "num_sub_spaces"):
-        target_b_split = target_b.subfunctions
-        source_b_split = source_b.subfunctions
     else:
-        target_b_split = [target_b]
-        source_b_split = [source_b]
-
-    # Apply adjoint transfer operator to each component
-    for i, (t_b, s_b) in enumerate(zip(target_b_split, source_b_split)):
-        if transfer_method == "interpolate":
-            s_b.interpolate(t_b, adjoint=True, **kwargs)
-        elif transfer_method == "project":
-            ksp = petsc4py.KSP().create()
-            ksp.setOperators(assemble_mass_matrix(t_b.function_space(), lumped=bounded))
-            mixed_mass = assemble_mixed_mass_matrix(Vt[i], Vs[i])
-            with t_b.dat.vec_ro as tb, s_b.dat.vec_wo as sb:
-                residual = tb.copy()
-                ksp.solveTranspose(tb, residual)
-                mixed_mass.mult(residual, sb)  # NOTE: already transposed above
+        if hasattr(Vt, "num_sub_spaces"):
+            for s, t in zip(source.subfunctions, target.subfunctions):
+                Vs = s.function_space()
+                Vt = t.function_space()
+                ksp = petsc4py.KSP().create()
+                ksp.setOperators(assemble_mass_matrix(Vs, lumped=bounded))
+                mixed_mass = assemble_mixed_mass_matrix(Vs, Vt)
+                with s.dat.vec_ro as vs, t.dat.vec_wo as vt:
+                    residual = vs.copy()
+                    ksp.solveTranspose(vs, residual)
+                    mixed_mass.mult(residual, vt)  # NOTE: already transposed above
         else:
-            raise ValueError(
-                f"Invalid transfer method: {transfer_method}."
-                " Options are 'interpolate' and 'project'."
-            )
-
-    # Map back to a Cofunction
-    return function2cofunction(source_b)
-
-
-def _validate_matching_spaces(Vs, Vt):
-    if hasattr(Vs, "num_sub_spaces"):
-        if not hasattr(Vt, "num_sub_spaces"):
-            raise ValueError(
-                "Source space has multiple components but target space does not."
-            )
-        if Vs.num_sub_spaces() != Vt.num_sub_spaces():
-            raise ValueError(
-                "Inconsistent numbers of components in source and target spaces:"
-                f" {Vs.num_sub_spaces()} vs. {Vt.num_sub_spaces()}."
-            )
-    elif hasattr(Vt, "num_sub_spaces"):
-        raise ValueError(
-            "Target space has multiple components but source space does not."
-        )
+            Vs = source.function_space()
+            Vt = target.function_space()
+            ksp = petsc4py.KSP().create()
+            ksp.setOperators(assemble_mass_matrix(Vs, lumped=bounded))
+            mixed_mass = assemble_mixed_mass_matrix(Vs, Vt)
+            with source.dat.vec_ro as vs, target.dat.vec_wo as vt:
+                residual = vs.copy()
+                ksp.solveTranspose(vs, residual)
+                mixed_mass.mult(residual, vt)  # NOTE: already transposed above
 
 
 @PETSc.Log.EventDecorator()

From 7b619028fe995632d12f57e87050ff2c573d9376 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:29:58 +0100
Subject: [PATCH 04/11] Subfunctions always available

---
 animate/interpolation.py | 31 +++++++------------------------
 1 file changed, 7 insertions(+), 24 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index 7b4475ce..b915a7eb 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -162,36 +162,19 @@ def project(source, target, **kwargs):
     """
     _validate_consistent_spaces(source.function_space(), target.function_space())
     if not source.function_space()._dual:
-        if hasattr(target.function_space(), "num_sub_spaces"):
-            for s, t in zip(source.subfunctions, target.subfunctions):
-                if bounded:
-                    _supermesh_project(s, t, bounded=True)
-                else:
-                    t.project(s, **kwargs)
-        else:
+        for s, t in zip(source.subfunctions, target.subfunctions):
             if bounded:
-                _supermesh_project(source, target, bounded=True)
+                _supermesh_project(s, t, bounded=True)
             else:
-                target.project(source, **kwargs)
+                t.project(s, **kwargs)
     else:
-        if hasattr(Vt, "num_sub_spaces"):
-            for s, t in zip(source.subfunctions, target.subfunctions):
-                Vs = s.function_space()
-                Vt = t.function_space()
-                ksp = petsc4py.KSP().create()
-                ksp.setOperators(assemble_mass_matrix(Vs, lumped=bounded))
-                mixed_mass = assemble_mixed_mass_matrix(Vs, Vt)
-                with s.dat.vec_ro as vs, t.dat.vec_wo as vt:
-                    residual = vs.copy()
-                    ksp.solveTranspose(vs, residual)
-                    mixed_mass.mult(residual, vt)  # NOTE: already transposed above
-        else:
-            Vs = source.function_space()
-            Vt = target.function_space()
+        for s, t in zip(source.subfunctions, target.subfunctions):
+            Vs = s.function_space()
+            Vt = t.function_space()
             ksp = petsc4py.KSP().create()
             ksp.setOperators(assemble_mass_matrix(Vs, lumped=bounded))
             mixed_mass = assemble_mixed_mass_matrix(Vs, Vt)
-            with source.dat.vec_ro as vs, target.dat.vec_wo as vt:
+            with s.dat.vec_ro as vs, t.dat.vec_wo as vt:
                 residual = vs.copy()
                 ksp.solveTranspose(vs, residual)
                 mixed_mass.mult(residual, vt)  # NOTE: already transposed above

From 21a2c9bf5c2d7540e751a2ee904992c8590b8b93 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:45:47 +0100
Subject: [PATCH 05/11] Add primal/dual consistency check

---
 animate/interpolation.py   | 6 ++----
 test/test_interpolation.py | 7 +++++++
 2 files changed, 9 insertions(+), 4 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index b915a7eb..acb34e49 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -73,6 +73,8 @@ def transfer(source, target_space, transfer_method="project", **kwargs):
 
 
 def _validate_consistent_spaces(Vs, Vt):
+    if Vs._dual != Vt._dual:
+        raise ValueError("Spaces must be both primal or both dual.")
     if hasattr(Vs, "num_sub_spaces"):
         if not hasattr(Vt, "num_sub_spaces"):
             raise ValueError(
@@ -87,10 +89,6 @@ def _validate_consistent_spaces(Vs, Vt):
         raise ValueError(
             "Target space has multiple components but source space does not."
         )
-    if Vs._dual and not Vt._dual:
-        raise ValueError("Source space is a dual space but target space is not.")
-    if not Vs._dual and Vt._dual:
-        raise ValueError("Target space is a dual space but source space is not.")
 
 
 @PETSc.Log.EventDecorator()
diff --git a/test/test_interpolation.py b/test/test_interpolation.py
index 485937c8..95e74ce1 100644
--- a/test/test_interpolation.py
+++ b/test/test_interpolation.py
@@ -186,6 +186,13 @@ def test_notimplemented_error(self, transfer_method):
         msg = f"Can only currently {transfer_method} Functions and Cofunctions."
         self.assertEqual(str(cm.exception), msg)
 
+    def test_primal_dual_inconsistency_error(self):
+        Vs = FunctionSpace(self.source_mesh, "CG", 1)
+        Vt = FunctionSpace(self.target_mesh, "CG", 1).dual()
+        with self.assertRaises(ValueError) as cm:
+            transfer(Function(Vs), Function(Vt))
+        self.assertEqual(str(cm.exception), "Spaces must be both primal or both dual.")
+
     @parameterized.expand(["interpolate", "project"])
     def test_no_sub_source_space(self, transfer_method):
         Vs = FunctionSpace(self.source_mesh, "CG", 1)

From 5cb4cd9738d4c5fbbdbdc027876f416c5e1c590d Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:46:13 +0100
Subject: [PATCH 06/11] Misc fixes and test updates

---
 animate/interpolation.py   |  2 +-
 test/test_interpolation.py | 22 ++++++++++------------
 2 files changed, 11 insertions(+), 13 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index acb34e49..300f5bff 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -138,7 +138,7 @@ def _supermesh_project(source, target, bounded=False):
 
 
 @PETSc.Log.EventDecorator()
-def project(source, target, **kwargs):
+def project(source, target, bounded=False, **kwargs):
     r"""
     Overload :func:`firedrake.projection.project` to account for the case of mixed
     function spaces.
diff --git a/test/test_interpolation.py b/test/test_interpolation.py
index 95e74ce1..0dcffa1c 100644
--- a/test/test_interpolation.py
+++ b/test/test_interpolation.py
@@ -22,8 +22,6 @@
 
 from animate.interpolation import (
     _supermesh_project,
-    _transfer_adjoint,
-    _transfer_forward,
     clement_interpolant,
     project,
     transfer,
@@ -168,12 +166,11 @@ def test_method_typo_error(self):
         msg = "Invalid transfer method: proj. Options are 'interpolate' and 'project'."
         self.assertEqual(str(cm.exception), msg)
 
-    @parameterized.expand([_transfer_forward, _transfer_adjoint])
-    def test_method_typo_error_private(self, private_transfer_func):
+    def test_method_typo_error(self):
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
         Vt = FunctionSpace(self.target_mesh, "CG", 1)
         with self.assertRaises(ValueError) as cm:
-            private_transfer_func(Function(Vs), Function(Vt), transfer_method="proj")
+            transfer(Function(Vs), Function(Vt), transfer_method="proj")
         msg = "Invalid transfer method: proj. Options are 'interpolate' and 'project'."
         self.assertEqual(str(cm.exception), msg)
 
@@ -210,7 +207,7 @@ def test_no_sub_source_space_adjoint(self, transfer_method):
         Vt = Vt * Vt
         with self.assertRaises(ValueError) as cm:
             transfer(Cofunction(Vs.dual()), Cofunction(Vt.dual()), transfer_method)
-        msg = "Source space has multiple components but target space does not."
+        msg = "Target space has multiple components but source space does not."
         self.assertEqual(str(cm.exception), msg)
 
     @parameterized.expand(["interpolate", "project"])
@@ -230,7 +227,7 @@ def test_no_sub_target_space_adjoint(self, transfer_method):
         Vs = Vs * Vs
         with self.assertRaises(ValueError) as cm:
             transfer(Cofunction(Vs.dual()), Cofunction(Vt.dual()), transfer_method)
-        msg = "Target space has multiple components but source space does not."
+        msg = "Source space has multiple components but target space does not."
         self.assertEqual(str(cm.exception), msg)
 
     @parameterized.expand(["interpolate", "project"])
@@ -384,11 +381,11 @@ def check_no_new_extrema(source, target, tol=1.0e-08):
     @parameterized.expand([(True, True), (True, False), (False, True), (False, False)])
     def test_supermesh_project(self, same_mesh, same_degree, tol=1.0e-07):
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
+        x, y = ufl.SpatialCoordinate(self.source_mesh)
+        source = Function(Vs).interpolate(self.sinusoid())
         target_mesh = self.source_mesh if same_mesh else self.target_mesh
         target_degree = 1 if same_degree else 0
         Vt = FunctionSpace(target_mesh, "DG", target_degree)
-        x, y = ufl.SpatialCoordinate(self.source_mesh)
-        source = Function(Vs).interpolate(self.sinusoid())
         target = Function(Vt)
         _supermesh_project(source, target, bounded=False)
         expected = Function(Vt).project(source)
@@ -400,10 +397,11 @@ def test_supermesh_project(self, same_mesh, same_degree, tol=1.0e-07):
     @parameterized.expand([(True,), (False,)])
     def test_mass_lumping(self, same_mesh, tol=1.0e-08):
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
-        target_mesh = self.source_mesh if same_mesh else self.target_mesh
-        Vt = FunctionSpace(target_mesh, "CG", 1)
         x, y = ufl.SpatialCoordinate(self.source_mesh)
         source = Function(Vs).interpolate(self.sinusoid())
-        target = project(source, Vt, bounded=True)
+        target_mesh = self.source_mesh if same_mesh else self.target_mesh
+        Vt = FunctionSpace(target_mesh, "CG", 1)
+        target = Function(Vt)
+        project(source, target, bounded=True)
         self.assertTrue(self.check_conservation(source, target, tol=tol))
         self.assertTrue(self.check_no_new_extrema(source, target, tol=tol))

From 98965716850f3e486f30a16bf7ce469f83b38f00 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:52:06 +0100
Subject: [PATCH 07/11] Adjoint interp tests for same spaces working

---
 test/test_interpolation.py | 24 ++++++++++--------------
 1 file changed, 10 insertions(+), 14 deletions(-)

diff --git a/test/test_interpolation.py b/test/test_interpolation.py
index 0dcffa1c..83e7b1db 100644
--- a/test/test_interpolation.py
+++ b/test/test_interpolation.py
@@ -26,7 +26,7 @@
     project,
     transfer,
 )
-from animate.utility import function2cofunction
+from animate.utility import cofunction2function, function2cofunction
 
 
 class TestClement(unittest.TestCase):
@@ -261,14 +261,12 @@ def test_transfer_same_space(self, transfer_method):
 
     @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_space_adjoint(self, transfer_method):
-        pytest.skip()  # TODO: (#114)
         Vs = FunctionSpace(self.source_mesh, "CG", 1)
-        source = Function(Vs).interpolate(self.sinusoid())
-        source = function2cofunction(source)
+        expected = Function(Vs).interpolate(self.sinusoid())
+        source = function2cofunction(expected)
         target = Cofunction(Vs.dual())
         transfer(source, target, transfer_method)
-        expected = source
-        self.assertAlmostEqual(errornorm(expected, target), 0)
+        self.assertAlmostEqual(errornorm(expected, cofunction2function(target)), 0)
 
     @parameterized.expand(["project", "interpolate"])
     def test_transfer_same_space_mixed(self, transfer_method):
@@ -285,18 +283,16 @@ def test_transfer_same_space_mixed(self, transfer_method):
 
     @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_space_mixed_adjoint(self, transfer_method):
-        pytest.skip()  # TODO: (#114)
         P1 = FunctionSpace(self.source_mesh, "CG", 1)
         Vs = P1 * P1
-        source = Function(Vs)
-        s1, s2 = source.subfunctions
-        s1.interpolate(self.sinusoid())
-        s2.interpolate(-self.sinusoid())
-        source = function2cofunction(source)
+        expected = Function(Vs)
+        e1, e2 = expected.subfunctions
+        e1.interpolate(self.sinusoid())
+        e2.interpolate(-self.sinusoid())
+        source = function2cofunction(expected)
         target = Cofunction(Vs.dual())
         transfer(source, target, transfer_method)
-        expected = source
-        self.assertAlmostEqual(errornorm(expected, target), 0)
+        self.assertAlmostEqual(errornorm(expected, cofunction2function(target)), 0)
 
     @parameterized.expand(["interpolate", "project"])
     def test_transfer_same_mesh(self, transfer_method):

From aa301f22348676996989164ed12e45253c1b719b Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 08:53:55 +0100
Subject: [PATCH 08/11] Attempt to fix project case

---
 animate/interpolation.py | 10 ++++++----
 1 file changed, 6 insertions(+), 4 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index 300f5bff..c1339203 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -167,15 +167,17 @@ def project(source, target, bounded=False, **kwargs):
                 t.project(s, **kwargs)
     else:
         for s, t in zip(source.subfunctions, target.subfunctions):
-            Vs = s.function_space()
-            Vt = t.function_space()
+            sf = cofunction2function(s)
+            tf = cofunction2function(t)
+            Vs = sf.function_space()
             ksp = petsc4py.KSP().create()
             ksp.setOperators(assemble_mass_matrix(Vs, lumped=bounded))
-            mixed_mass = assemble_mixed_mass_matrix(Vs, Vt)
-            with s.dat.vec_ro as vs, t.dat.vec_wo as vt:
+            mixed_mass = assemble_mixed_mass_matrix(Vs, tf.function_space())
+            with sf.dat.vec_ro as vs, tf.dat.vec_wo as vt:
                 residual = vs.copy()
                 ksp.solveTranspose(vs, residual)
                 mixed_mass.mult(residual, vt)  # NOTE: already transposed above
+            function2cofunction(tf, cofunc=t)
 
 
 @PETSc.Log.EventDecorator()

From 9f8b6bbee0b611b809e53800106a22f220c50291 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Thu, 29 May 2025 09:04:22 +0100
Subject: [PATCH 09/11] Lint

---
 test/test_interpolation.py | 8 --------
 1 file changed, 8 deletions(-)

diff --git a/test/test_interpolation.py b/test/test_interpolation.py
index 83e7b1db..f79b6c30 100644
--- a/test/test_interpolation.py
+++ b/test/test_interpolation.py
@@ -166,14 +166,6 @@ def test_method_typo_error(self):
         msg = "Invalid transfer method: proj. Options are 'interpolate' and 'project'."
         self.assertEqual(str(cm.exception), msg)
 
-    def test_method_typo_error(self):
-        Vs = FunctionSpace(self.source_mesh, "CG", 1)
-        Vt = FunctionSpace(self.target_mesh, "CG", 1)
-        with self.assertRaises(ValueError) as cm:
-            transfer(Function(Vs), Function(Vt), transfer_method="proj")
-        msg = "Invalid transfer method: proj. Options are 'interpolate' and 'project'."
-        self.assertEqual(str(cm.exception), msg)
-
     @parameterized.expand(["interpolate", "project"])
     def test_notimplemented_error(self, transfer_method):
         Vs = FunctionSpace(self.source_mesh, "CG", 1)

From 6de0c0a40083f654b0ffd5801199cadb3a29009c Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Mon, 17 Nov 2025 07:19:35 +0000
Subject: [PATCH 10/11] Replace __future__ with interpolation

---
 animate/interpolation.py | 10 +++++-----
 animate/metric.py        |  2 +-
 animate/quality.py       |  2 +-
 animate/recovery.py      |  2 +-
 animate/utility.py       |  2 +-
 5 files changed, 9 insertions(+), 9 deletions(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index d5578efa..440a08e9 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -24,7 +24,7 @@
 @PETSc.Log.EventDecorator()
 def transfer(source, target_space, transfer_method="project", **kwargs):
     r"""
-    Overload functions :func:`firedrake.__future__.interpolate` and
+    Overload functions :func:`firedrake.interpolation.interpolate` and
     :func:`firedrake.projection.project` to account for the case of two mixed
     function spaces defined on different meshes and for the adjoint interpolation
     operator when applied to :class:`firedrake.cofunction.Cofunction`\s.
@@ -46,7 +46,7 @@ def transfer(source, target_space, transfer_method="project", **kwargs):
     :rtype: :class:`firedrake.function.Function` or
         :class:`firedrake.cofunction.Cofunction`
 
-    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate` or
+    Extra keyword arguments are passed to :func:`firedrake.interpolation.interpolate` or
         :func:`firedrake.projection.project`.
     """
     if transfer_method not in ("interpolate", "project"):
@@ -94,8 +94,8 @@ def _validate_consistent_spaces(Vs, Vt):
 @PETSc.Log.EventDecorator()
 def interpolate(source, target, **kwargs):
     r"""
-    Overload :func:`firedrake.__future__.interpolate` to account for the case of mixed
-    function spaces.
+    Overload :func:`firedrake.interpolation.interpolate` to account for the case of
+    mixed function spaces.
 
     :arg source: the function or cofunction to be transferred
     :type source: :class:`firedrake.function.Function` or
@@ -105,7 +105,7 @@ def interpolate(source, target, **kwargs):
     :type target: :class:`firedrake.function.Function` or
         :class:`firedrake.cofunction.Cofunction`
 
-    Extra keyword arguments are passed to :func:`firedrake.__future__.interpolate`
+    Extra keyword arguments are passed to :func:`firedrake.interpolation.interpolate`
     """
     _validate_consistent_spaces(source.function_space(), target.function_space())
     if hasattr(target.function_space(), "num_sub_spaces"):
diff --git a/animate/metric.py b/animate/metric.py
index b6702d32..2c2be202 100644
--- a/animate/metric.py
+++ b/animate/metric.py
@@ -7,7 +7,7 @@
 import numpy as np
 import sympy
 import ufl
-from firedrake.__future__ import interpolate
+from firedrake.interpolation import interpolate
 from firedrake.petsc import PETSc
 from pyop2 import op2
 
diff --git a/animate/quality.py b/animate/quality.py
index d3e060f4..fc47226c 100644
--- a/animate/quality.py
+++ b/animate/quality.py
@@ -6,7 +6,7 @@
 
 import firedrake
 import ufl
-from firedrake.__future__ import interpolate
+from firedrake.interpolation import interpolate
 from firedrake.petsc import PETSc
 from pyop2 import op2
 from pyop2.utils import get_petsc_dir
diff --git a/animate/recovery.py b/animate/recovery.py
index eb833589..051c9507 100644
--- a/animate/recovery.py
+++ b/animate/recovery.py
@@ -6,7 +6,7 @@
 
 import firedrake
 import ufl
-from firedrake.__future__ import interpolate
+from firedrake.interpolation import interpolate
 from firedrake.petsc import PETSc
 from pyop2 import op2
 
diff --git a/animate/utility.py b/animate/utility.py
index 3f9ad5cd..0dacf50a 100644
--- a/animate/utility.py
+++ b/animate/utility.py
@@ -8,7 +8,7 @@
 import firedrake.function as ffunc
 import firedrake.mesh as fmesh
 import ufl
-from firedrake.__future__ import interpolate
+from firedrake.interpolation import interpolate
 from firedrake.petsc import PETSc
 from mpi4py import MPI
 

From fe602c3c6115952a341f452cc0327ddda8056989 Mon Sep 17 00:00:00 2001
From: Joe Wallwork <jw2423@cam.ac.uk>
Date: Mon, 17 Nov 2025 07:21:03 +0000
Subject: [PATCH 11/11] Mention dual functions and boundedness in project
 docstring

---
 animate/interpolation.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/animate/interpolation.py b/animate/interpolation.py
index 440a08e9..5f77102f 100644
--- a/animate/interpolation.py
+++ b/animate/interpolation.py
@@ -141,7 +141,7 @@ def _supermesh_project(source, target, bounded=False):
 def project(source, target, bounded=False, **kwargs):
     r"""
     Overload :func:`firedrake.projection.project` to account for the case of mixed
-    function spaces.
+    function spaces and handle dual Functions, with a bounded option.
 
     For details on the approach for achieving boundedness through mass lumping and
     post-processing, see :cite:`Farrell:2009`.