FEniCS / dolfinx

Next generation FEniCS problem solving environment
https://fenicsproject.org
GNU Lesser General Public License v3.0
698 stars 172 forks source link

Revise layout of coordinate data in Python interface #2541

Open garth-wells opened 1 year ago

garth-wells commented 1 year ago

We have a bit of a muddle with the layout of coordinate data in the Python interface. The background is:

  1. In the C++ layer we normally work with shape (num_points, gdim) and row-major storage
  2. From Python it is often natural to have shape (gdim, num_points). The allow a natural syntax for vectorised NumPy operations by coordinate component, e.g. x[1] *= 2.0 multiplies the y-component of all points by 2.

Should we make the layout consistent in the Python interface? Should the shape be (gdim, num_points) throughout?

See related issue #2535.

jorgensd commented 1 year ago

I guess the issue is that a lot of the functions in C++ are not vectorized. It would mean a lot of strided views of data in: https://github.com/FEniCS/dolfinx/blob/main/cpp/dolfinx/geometry/utils.h https://github.com/FEniCS/dolfinx/blob/main/cpp/dolfinx/geometry/gjk.h I guess it will have an effect for assembly as well, however, there we already do a strided copy: https://github.com/FEniCS/dolfinx/blob/main/cpp/dolfinx/fem/assemble_scalar_impl.h#L49-L54

garth-wells commented 1 year ago

I guess the issue is that a lot of the functions in C++ are not vectorized.

The main change would be copying of data between layouts inside the pybind11 layer, but not necessarily in performance critical parts.