Fix for Cofunction self-assignment via interpolation

firedrakeproject · Dec 18, 2024 · b8c5252 · b8c5252
1 parent 18a952f
commit b8c5252
Show file tree

Hide file tree

Showing 2 changed files with 13 additions and 1 deletion.
diff --git a/firedrake/interpolation.py b/firedrake/interpolation.py
@@ -859,7 +859,12 @@ def _interpolate(self, *function, output=None, transpose=False, **kwargs):
                 V = self.V
             result = output or firedrake.Function(V)
             with function.dat.vec_ro as x, result.dat.vec_wo as out:
-                mul(x, out)
+                if x.id != out.id:
+                    mul(x, out)
+                else:
+                    out_ = out.duplicate()
+                    mul(x, out_)
+                    out_.copy(result=out)
             return result
 
         else:

diff --git a/tests/firedrake/regression/test_interp_dual.py b/tests/firedrake/regression/test_interp_dual.py
@@ -34,6 +34,13 @@ def f1(mesh, V1):
     return Function(V1).interpolate(expr)
 
 
+def test_interp_self(V1):
+    a = assemble(TestFunction(V1) * dx)
+    b = assemble(TestFunction(V1) * dx)
+    a.interpolate(a)
+    assert (a.dat.data_ro == b.dat.data_ro).all()
+
+
 def test_assemble_interp_operator(V2, f1):
     # Check type
     If1 = Interpolate(f1, V2)