Add void* to tabulate tensor kernel

.github/workflows/ccpp.yml

@@ -115,13 +115,13 @@ jobs:
         run: |
           pip install git+https://github.com/fenics/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/fenics/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/fenics/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
       - name: Install FEniCS Python components
         if: github.event_name == 'workflow_dispatch'
         run: |
           pip install git+https://github.com/FEniCS/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/FEniCS/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/FEniCS/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
 
       - name: Configure and install C++
         run: |
@@ -188,14 +188,14 @@ jobs:
         run: |
           pip install git+https://github.com/FEniCS/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/FEniCS/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/FEniCS/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
 
       - name: Install FEniCS Python components
         if: github.event_name == 'workflow_dispatch'
         run: |
           pip install git+https://github.com/FEniCS/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/FEniCS/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/FEniCS/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
 
       - name: Configure, build and install C++ library
         run: |
@@ -270,13 +270,13 @@ jobs:
         run: |
           pip install git+https://github.com/FEniCS/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/FEniCS/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/FEniCS/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
       - name: Install FEniCS Python components
         if: github.event_name == 'workflow_dispatch'
         run: |
           pip install git+https://github.com/FEniCS/ufl.git@${{ env.ufl_ref }}
           pip install git+https://github.com/FEniCS/basix.git@${{ env.basix_ref }}
-          pip install git+https://github.com/FEniCS/ffcx.git@${{ env.ffcx_ref }}
+          pip install git+https://github.com/sclaus2/ffcx.git@${{ env.ffcx_ref }}
 
       - name: Configure C++
         run: |

.github/workflows/fenicsx-refs.env

@@ -1,3 +1,3 @@
 basix_ref=main
 ufl_ref=main
-ffcx_ref=main
+ffcx_ref=sclaus/add-void-to-kernels

cpp/demo/custom_kernel/main.cpp

@@ -238,7 +238,7 @@ void assemble(MPI_Comm comm)
   std::array<T, 9> A_hat_b = A_ref<T>(phi, weights);
   auto kernel_a
       = [A_hat = mdspan2_t<T, 3, 3>(A_hat_b.data()),
-         detJ](T* A, const T*, const T*, const T* x, const int*, const uint8_t*)
+         detJ](T* A, const T*, const T*, const T* x, const int*, const uint8_t*, void*)
   {
     T scale = detJ(mdspan2_t<const T, 3, 3>(x));
     mdspan2_t<T, 3, 3> _A(A);
@@ -250,7 +250,7 @@ void assemble(MPI_Comm comm)
   // Finite element RHS (f=1) kernel function
   auto kernel_L
       = [b_hat = b_ref<T>(phi, weights),
-         detJ](T* b, const T*, const T*, const T* x, const int*, const uint8_t*)
+         detJ](T* b, const T*, const T*, const T* x, const int*, const uint8_t*, void*)
   {
     T scale = detJ(mdspan2_t<const T, 3, 3>(x));
     for (std::size_t i = 0; i < 3; ++i)

cpp/dolfinx/fem/Expression.h

@@ -74,7 +74,7 @@ class Expression
           constants,
       std::span<const geometry_type> X, std::array<std::size_t, 2> Xshape,
       std::function<void(scalar_type*, const scalar_type*, const scalar_type*,
-                         const geometry_type*, const int*, const uint8_t*)>
+                         const geometry_type*, const int*, const uint8_t*, void*)>
           fn,
       const std::vector<std::size_t>& value_shape,
       std::shared_ptr<const FunctionSpace<geometry_type>> argument_space
@@ -144,7 +144,7 @@ class Expression
 
   /// @brief Function for tabulating the Expression.
   const std::function<void(scalar_type*, const scalar_type*, const scalar_type*,
-                           const geometry_type*, const int*, const uint8_t*)>&
+                           const geometry_type*, const int*, const uint8_t*, void*)>&
   kernel() const
   {
     return _fn;
@@ -179,7 +179,7 @@ class Expression
 
   // Function to evaluate the Expression
   std::function<void(scalar_type*, const scalar_type*, const scalar_type*,
-                     const geometry_type*, const int*, const uint8_t*)>
+                     const geometry_type*, const int*, const uint8_t*, void*)>
       _fn;
 
   // Shape of the evaluated Expression

cpp/dolfinx/fem/Form.h

@@ -117,7 +117,7 @@ struct integral_data
     requires std::is_convertible_v<
                  std::remove_cvref_t<K>,
                  std::function<void(T*, const T*, const T*, const U*,
-                                    const int*, const uint8_t*)>>
+                                    const int*, const uint8_t*, void*)>>
                  and std::is_convertible_v<std::remove_cvref_t<V>,
                                            std::vector<std::int32_t>>
                  and std::is_convertible_v<std::remove_cvref_t<W>,
@@ -130,7 +130,7 @@ struct integral_data
 
   /// @brief The integration kernel.
   std::function<void(T*, const T*, const T*, const U*, const int*,
-                     const uint8_t*)>
+                     const uint8_t*, void*)>
       kernel;
 
   /// @brief The entities to integrate over for this integral. These are
@@ -390,7 +390,7 @@ class Form
   /// kernels for a given ID in mixed-topology meshes).
   /// @return Function to call for `tabulate_tensor`.
   std::function<void(scalar_type*, const scalar_type*, const scalar_type*,
-                     const geometry_type*, const int*, const uint8_t*)>
+                     const geometry_type*, const int*, const uint8_t*, void*)>
   kernel(IntegralType type, int id, int kernel_idx) const
   {
     auto it = _integrals.find({type, id, kernel_idx});

cpp/dolfinx/fem/assemble_expression_impl.h

@@ -92,7 +92,7 @@ void tabulate_expression(
                     std::next(coord_dofs.begin(), 3 * i));
       }
       fn(values_local.data(), &coeffs(e, 0), constants.data(),
-         coord_dofs.data(), nullptr, nullptr);
+         coord_dofs.data(), nullptr, nullptr, nullptr);
     }
     else
     {
@@ -104,7 +104,7 @@ void tabulate_expression(
                     std::next(coord_dofs.begin(), 3 * i));
       }
       fn(values_local.data(), &coeffs(e, 0), constants.data(),
-         coord_dofs.data(), &entities(e, 1), nullptr);
+         coord_dofs.data(), &entities(e, 1), nullptr, nullptr);
     }
 
     P0(values_local, cell_info, e, size0);

cpp/dolfinx/fem/assemble_matrix_impl.h

@@ -110,7 +110,7 @@ void assemble_cells(
     // Tabulate tensor
     std::ranges::fill(Ae, 0);
     kernel(Ae.data(), &coeffs(c, 0), constants.data(), cdofs.data(), nullptr,
-           nullptr);
+           nullptr, nullptr);
 
     // Compute A = P_0 \tilde{A} P_1^T (dof transformation)
     P0(_Ae, cell_info0, cell0, ndim1);  // B = P0 \tilde{A}
@@ -253,7 +253,7 @@ void assemble_exterior_facets(
     // Tabulate tensor
     std::ranges::fill(Ae, 0);
     kernel(Ae.data(), &coeffs(f, 0), constants.data(), cdofs.data(),
-           &local_facet, &perm);
+           &local_facet, &perm, nullptr);
 
     P0(_Ae, cell_info0, cell0, ndim1);
     P1T(_Ae, cell_info1, cell1, ndim0);
@@ -423,7 +423,7 @@ void assemble_interior_facets(
                           : std::array{perms(cells[0], local_facet[0]),
                                        perms(cells[1], local_facet[1])};
     kernel(Ae.data(), &coeffs(f, 0, 0), constants.data(), cdofs.data(),
-           local_facet.data(), perm.data());
+           local_facet.data(), perm.data(), nullptr);
 
     // Local element layout is a 2x2 block matrix with structure
     //

cpp/dolfinx/fem/assemble_scalar_impl.h

@@ -51,7 +51,7 @@ T assemble_cells(mdspan2_t x_dofmap,
     }
 
     fn(&value, &coeffs(index, 0), constants.data(), cdofs.data(), nullptr,
-       nullptr);
+       nullptr, nullptr);
   }
 
   return value;
@@ -92,7 +92,7 @@ T assemble_exterior_facets(
     // Permutations
     std::uint8_t perm = perms.empty() ? 0 : perms(cell, local_facet);
     fn(&value, &coeffs(f, 0), constants.data(), cdofs.data(), &local_facet,
-       &perm);
+       &perm, nullptr);
   }
 
   return value;
@@ -144,7 +144,7 @@ T assemble_interior_facets(
                           : std::array{perms(cells[0], local_facet[0]),
                                        perms(cells[1], local_facet[1])};
     fn(&value, &coeffs(f, 0, 0), constants.data(), cdofs.data(),
-       local_facet.data(), perm.data());
+       local_facet.data(), perm.data(), nullptr);
   }
 
   return value;

cpp/dolfinx/fem/assemble_vector_impl.h

@@ -163,7 +163,7 @@ void _lift_bc_cells(
     Ae.resize(num_rows * num_cols);
     std::ranges::fill(Ae, 0);
     kernel(Ae.data(), &coeffs(index, 0), constants.data(), cdofs.data(),
-           nullptr, nullptr);
+           nullptr, nullptr, nullptr);
     P0(Ae, cell_info0, c0, num_cols);
     P1T(Ae, cell_info1, c1, num_rows);
 
@@ -345,7 +345,7 @@ void _lift_bc_exterior_facets(
     Ae.resize(num_rows * num_cols);
     std::ranges::fill(Ae, 0);
     kernel(Ae.data(), &coeffs(index, 0), constants.data(), cdofs.data(),
-           &local_facet, &perm);
+           &local_facet, &perm, nullptr);
     P0(Ae, cell_info0, cell0, num_cols);
     P1T(Ae, cell_info1, cell1, num_rows);
 
@@ -544,7 +544,7 @@ void _lift_bc_interior_facets(
                           : std::array{perms(cells[0], local_facet[0]),
                                        perms(cells[1], local_facet[1])};
     kernel(Ae.data(), &coeffs(f, 0, 0), constants.data(), cdofs.data(),
-           local_facet.data(), perm.data());
+           local_facet.data(), perm.data(), nullptr);
 
     std::span<T> _Ae(Ae);
     std::span<T> sub_Ae0 = _Ae.subspan(bs0 * dmap0_cell0.size() * num_cols,
@@ -676,7 +676,7 @@ void assemble_cells(
     // Tabulate vector for cell
     std::ranges::fill(be, 0);
     kernel(be.data(), &coeffs(index, 0), constants.data(), cdofs.data(),
-           nullptr, nullptr);
+           nullptr, nullptr, nullptr);
     P0(_be, cell_info0, c0, 1);
 
     // Scatter cell vector to 'global' vector array
@@ -769,7 +769,7 @@ void assemble_exterior_facets(
     // Tabulate element vector
     std::ranges::fill(be, 0);
     fn(be.data(), &coeffs(f, 0), constants.data(), cdofs.data(), &local_facet,
-       &perm);
+       &perm, nullptr);
 
     P0(_be, cell_info0, cell0, 1);
 
@@ -877,7 +877,7 @@ void assemble_interior_facets(
                           : std::array{perms(cells[0], local_facet[0]),
                                        perms(cells[1], local_facet[1])};
     fn(be.data(), &coeffs(f, 0, 0), constants.data(), cdofs.data(),
-       local_facet.data(), perm.data());
+       local_facet.data(), perm.data(), nullptr);
 
     std::span<T> _be(be);
     std::span<T> sub_be = _be.subspan(bs * dmap0.size(), bs * dmap1.size());

cpp/dolfinx/fem/traits.h

@@ -28,7 +28,7 @@ concept DofTransformKernel
 template <class U, class T>
 concept FEkernel = std::is_invocable_v<U, T*, const T*, const T*,
                                        const scalar_value_t<T>*,
-                                       const int*, const std::uint8_t*>;
+                                       const int*, const std::uint8_t*, void*>;
 
 /// @brief Concept for mdspan of rank 1 or 2.
 template <class T>

cpp/dolfinx/fem/utils.h

@@ -447,7 +447,7 @@ Form<T, U> create_form_factory(
   // Get list of integral IDs, and load tabulate tensor into memory for
   // each
   using kern_t = std::function<void(T*, const T*, const T*, const U*,
-                                    const int*, const std::uint8_t*)>;
+                                    const int*, const std::uint8_t*, void*)>;
   std::map<std::tuple<IntegralType, int, int>, integral_data<T, U>> integrals;
 
   auto check_geometry_hash
@@ -498,7 +498,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex64);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex64);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
         else if constexpr (std::is_same_v<T, double>)
@@ -508,7 +508,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex128);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex128);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
 
@@ -583,7 +583,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex64);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex64);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
         else if constexpr (std::is_same_v<T, double>)
@@ -593,7 +593,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex128);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex128);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
         assert(k);
@@ -689,7 +689,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex64);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex64);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
         else if constexpr (std::is_same_v<T, double>)
@@ -699,7 +699,7 @@ Form<T, U> create_form_factory(
         {
           k = reinterpret_cast<void (*)(
               T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-              const unsigned char*)>(integral->tabulate_tensor_complex128);
+              const unsigned char*, void*)>(integral->tabulate_tensor_complex128);
         }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
         assert(k);
@@ -918,7 +918,7 @@ Expression<T, U> create_expression(
   std::vector<std::size_t> value_shape(e.value_shape,
                                        e.value_shape + e.num_components);
   std::function<void(T*, const T*, const T*, const scalar_value_t<T>*,
-                     const int*, const std::uint8_t*)>
+                     const int*, const std::uint8_t*, void*)>
       tabulate_tensor = nullptr;
   if constexpr (std::is_same_v<T, float>)
     tabulate_tensor = e.tabulate_tensor_float32;
@@ -927,7 +927,7 @@ Expression<T, U> create_expression(
   {
     tabulate_tensor = reinterpret_cast<void (*)(
         T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-        const unsigned char*)>(e.tabulate_tensor_complex64);
+        const unsigned char*, void*)>(e.tabulate_tensor_complex64);
   }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
   else if constexpr (std::is_same_v<T, double>)
@@ -937,7 +937,7 @@ Expression<T, U> create_expression(
   {
     tabulate_tensor = reinterpret_cast<void (*)(
         T*, const T*, const T*, const scalar_value_t<T>*, const int*,
-        const unsigned char*)>(e.tabulate_tensor_complex128);
+        const unsigned char*, void*)>(e.tabulate_tensor_complex128);
   }
 #endif // DOLFINX_NO_STDC_COMPLEX_KERNELS
   else