parameters = (Mesh(VectorElement(FiniteElement('Q', quadrilateral, 1), dim=3), 1670), Mesh(VectorElement(FiniteElement('Lagrange', ...FixedIndex(1),))))), array([ 1.        ,  1.        , -1.        , -1.        ,  1.41421356,
       -0.3660254 ]), ...)

    @pytest.mark.parallel
    def test_interpolate_cross_mesh_parallel(parameters):
>       test_interpolate_cross_mesh(parameters)

tests/firedrake/regression/test_interpolate_cross_mesh.py:741: 
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 
tests/firedrake/regression/test_interpolate_cross_mesh.py:302: in test_interpolate_cross_mesh
    get_expected_values(
tests/firedrake/regression/test_interpolate_cross_mesh.py:501: in get_expected_values
    interpolate_function(
tests/firedrake/regression/test_interpolate_cross_mesh.py:533: in interpolate_function
    got = f_dest.at(coords)
petsc4py/PETSc/Log.pyx:188: in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
    ???
petsc4py/PETSc/Log.pyx:189: in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
    ???
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 

self = Coefficient(WithGeometry(FunctionSpace(<firedrake.mesh.MeshTopology object at 0x7f682a6e3aa0>, FiniteElement('Disconti...s Lagrange', triangle, 2), name=None), Mesh(VectorElement(FiniteElement('Lagrange', triangle, 1), dim=3), 1671)), 3438)
arg = array([[ 0.        ,  1.        ,  0.        ],
       [ 1.        ,  0.        ,  0.        ],
       [-1.        ,  ...       ,  0.        ],
       [ 0.70710678,  0.70710678,  0.        ],
       [-0.8660254 ,  0.5       ,  0.        ]])
args = (), kwargs = {}
MPI = <module 'mpi4py.MPI' from '/__w/firedrake/firedrake_venv/lib/python3.12/site-packages/mpi4py/MPI.cpython-312-x86_64-linux-gnu.so'>
dont_raise = False
mesh = Mesh(VectorElement(FiniteElement('Lagrange', triangle, 1), dim=3), 1671)
gdim = 3
root_arg = array([[ 0.        ,  1.        ,  0.        ],
       [ 1.        ,  0.        ,  0.        ],
       [-1.        ,  ...       ,  0.        ],
       [ 0.70710678,  0.70710678,  0.        ],
       [-0.8660254 ,  0.5       ,  0.        ]])
same_arg = True, diff_arg = 0
single_eval = <function Function.at.<locals>.single_eval at 0x7f682bd65e40>
points = array([[ 0.        ,  1.        ,  0.        ],
       [ 1.        ,  0.        ,  0.        ],
       [-1.        ,  ...       ,  0.        ],
       [ 0.70710678,  0.70710678,  0.        ],
       [-0.8660254 ,  0.5       ,  0.        ]])
value_shape = ()

    @PETSc.Log.EventDecorator()
    def at(self, arg, *args, **kwargs):
        r"""Evaluate function at points.
    
        :arg arg: The point to locate.
        :arg args: Additional points.
        :kwarg dont_raise: Do not raise an error if a point is not found.
        :kwarg tolerance: Tolerence to use when checking if a point is
            in a cell. Default is the ``tolerance`` provided when
            creating the :func:`~.Mesh` the function is defined on.
            Changing this from default will cause the spatial index to
            be rebuilt which can take some time.
        """
        # Shortcut if function space is the R-space
        if self.ufl_element().family() == "Real":
            return self.dat.data_ro
    
        # Need to ensure data is up-to-date for reading
        self.dat.global_to_local_begin(op2.READ)
        self.dat.global_to_local_end(op2.READ)
        from mpi4py import MPI
    
        if args:
            arg = (arg,) + args
        arg = np.asarray(arg, dtype=utils.ScalarType)
        if utils.complex_mode:
            if not np.allclose(arg.imag, 0):
                raise ValueError("Provided points have non-zero imaginary part")
            arg = arg.real.copy()
    
        dont_raise = kwargs.get('dont_raise', False)
    
        tolerance = kwargs.get('tolerance', None)
        mesh = self.function_space().mesh()
        if tolerance is None:
            tolerance = mesh.tolerance
        else:
            mesh.tolerance = tolerance
    
        # Handle f.at(0.3)
        if not arg.shape:
            arg = arg.reshape(-1)
    
        if mesh.variable_layers:
            raise NotImplementedError("Point evaluation not implemented for variable layers")
    
        # Validate geometric dimension
        gdim = mesh.geometric_dimension()
        if arg.shape[-1] == gdim:
            pass
        elif len(arg.shape) == 1 and gdim == 1:
            arg = arg.reshape(-1, 1)
        else:
            raise ValueError("Point dimension (%d) does not match geometric dimension (%d)." % (arg.shape[-1], gdim))
    
        # Check if we have got the same points on each process
        root_arg = self._comm.bcast(arg, root=0)
        same_arg = arg.shape == root_arg.shape and np.allclose(arg, root_arg)
        diff_arg = self._comm.allreduce(int(not same_arg), op=MPI.SUM)
        if diff_arg:
            raise ValueError("Points to evaluate are inconsistent among processes.")
    
        def single_eval(x, buf):
            r"""Helper function to evaluate at a single point."""
            err = self._c_evaluate(tolerance=tolerance)(self._ctypes,
                                                        x.ctypes.data_as(POINTER(c_double)),
                                                        buf.ctypes.data_as(c_void_p))
            if err == -1:
                raise PointNotInDomainError(self.function_space().mesh(), x.reshape(-1))
    
        if not len(arg.shape) <= 2:
            raise ValueError("Function.at expects point or array of points.")
        points = arg.reshape(-1, arg.shape[-1])
        value_shape = self.ufl_shape
    
        subfunctions = self.subfunctions
        mixed = type(self.function_space().ufl_element()) is MixedElement
    
        # Local evaluation
        l_result = []
        for i, p in enumerate(points):
            try:
                if mixed:
                    l_result.append((i, tuple(f.at(p) for f in subfunctions)))
                else:
                    p_result = np.zeros(value_shape, dtype=ScalarType)
                    single_eval(points[i:i+1], p_result)
                    l_result.append((i, p_result))
            except PointNotInDomainError:
                # Skip point
                pass
    
        # Collecting the results
        def same_result(a, b):
            if mixed:
                for a_, b_ in zip(a, b):
                    if not np.allclose(a_, b_):
                        return False
                return True
            else:
                return np.allclose(a, b)
    
        all_results = self.comm.allgather(l_result)
        g_result = [None] * len(points)
        for results in all_results:
            for i, result in results:
                if g_result[i] is None:
                    g_result[i] = result
                elif same_result(result, g_result[i]):
                    pass
                else:
>                   raise RuntimeError("Point evaluation gave different results across processes.")
E                   RuntimeError: Point evaluation gave different results across processes.

firedrake/function.py:705: RuntimeError