From b8c5252ab3bb02335ea3939c59f0e35e160fc610 Mon Sep 17 00:00:00 2001 From: "James R. Maddison" Date: Wed, 18 Dec 2024 13:13:27 +0000 Subject: [PATCH] Fix for Cofunction self-assignment via interpolation --- firedrake/interpolation.py | 7 ++++++- tests/firedrake/regression/test_interp_dual.py | 7 +++++++ 2 files changed, 13 insertions(+), 1 deletion(-) diff --git a/firedrake/interpolation.py b/firedrake/interpolation.py index 0a598ba34b..449c4ff4d4 100644 --- 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 index 5928d81e4e..ac24e3239f 100644 --- 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)