Reverse indices for dual_arg

Author: pbrubeck

firedrake/assemble.py

@@ -596,14 +596,13 @@ def base_form_assembly_visitor(self, expr, tensor, *args):
             same_mesh = source_mesh.topology is target_mesh.topology
 
             # Get the interpolator
-            interp_data = expr.interp_data
+            interp_data = expr.interp_data.copy()
             default_missing_val = interp_data.pop('default_missing_val', None)
             if same_mesh and ((is_adjoint and rank == 1) or rank == 0):
-                interp_data = interp_data.copy()
                 interp_data["access"] = op2.INC
 
             dual_arg = v if same_mesh else V
-            interp_expr = firedrake.Interpolate(operand, v=dual_arg, **interp_data)
+            interp_expr = reconstruct_interp(operand, v=dual_arg)
             interpolator = firedrake.Interpolator(interp_expr, V, **interp_data)
 
             # Assembly

firedrake/interpolation.py

@@ -938,7 +938,7 @@ def callable():
 
             # Interpolate each sub expression into each function space
             for Vsub, sub_tensor, sub_op, sub_dual in zip(V, tensor, operands, duals):
-                sub_expr = ufl.Interpolate(sub_op, sub_dual)
+                sub_expr = expr._ufl_expr_reconstruct_(sub_op, sub_dual)
                 loops.extend(_interpolator(Vsub, sub_tensor, sub_expr, subset, arguments, access, bcs=bcs))
 
         if bcs and rank == 1:
@@ -1050,7 +1050,6 @@ def _interpolator(V, tensor, expr, subset, arguments, access, bcs=None):
     name = kernel.name
     kernel = op2.Kernel(ast, name, requires_zeroed_output_arguments=True,
                         flop_count=kernel.flop_count, events=(kernel.event,))
-
     parloop_args = [kernel, cell_set]
 
     expr = ufl.Interpolate(operand, v=dual_arg)
@@ -1069,6 +1068,7 @@ def _interpolator(V, tensor, expr, subset, arguments, access, bcs=None):
     else:
         copyin = ()
         copyout = ()
+
     if isinstance(tensor, op2.Global):
         parloop_args.append(tensor(access))
     elif isinstance(tensor, op2.Dat):

tsfc/driver.py

@@ -300,7 +300,7 @@ def compile_expression_dual_evaluation(expression, to_element, ufl_element, *,
 
     # Compute the adjoint by contracting against the dual argument
     if dual_arg in coefficients:
-        beta = basis_indices
+        beta = basis_indices[::-1]
         shape = tuple(i.extent for i in beta)
         gem_dual = gem.Variable(f"w_{coefficients.index(dual_arg)}", shape)
         evaluation = gem.IndexSum(evaluation * gem_dual[beta], beta)