From 1793cc36ab53dd6c121d4c4d0330b69f7a76a8ea Mon Sep 17 00:00:00 2001 From: "Garth N. Wells" Date: Sun, 9 Jun 2024 22:41:54 +0100 Subject: [PATCH] Replace typedefs --- cpp/Contact.cpp | 13 +++--- cpp/KernelData.cpp | 12 +++--- cpp/KernelData.h | 21 +++++----- cpp/QuadratureRule.cpp | 9 +++-- cpp/QuadratureRule.h | 7 +--- cpp/RayTracing.h | 4 +- cpp/coefficients.cpp | 12 +++--- cpp/coefficients.h | 6 +-- cpp/contact_kernels.cpp | 18 ++++++--- cpp/elasticity.cpp | 10 ++--- cpp/elasticity.h | 9 +++-- cpp/geometric_quantities.cpp | 10 +++-- cpp/geometric_quantities.h | 18 +++++---- cpp/meshtie_kernels.cpp | 15 ++++--- cpp/rigid_surface_kernels.cpp | 8 ++-- cpp/utils.cpp | 63 ++++++++++++++++------------- cpp/utils.h | 27 ++++++++----- python/dolfinx_contact/wrappers.cpp | 3 +- 18 files changed, 148 insertions(+), 117 deletions(-) diff --git a/cpp/Contact.cpp b/cpp/Contact.cpp index 43243101..30d79a1d 100644 --- a/cpp/Contact.cpp +++ b/cpp/Contact.cpp @@ -361,11 +361,12 @@ dolfinx_contact::Contact::pack_nx(int pair) const = dolfinx::mesh::cell_type_to_basix_type(topology->cell_type()); auto [facet_normalsb, n_shape] = basix::cell::facet_outward_normals(cell_type); - cmdspan2_t facet_normals(facet_normalsb.data(), n_shape); + mdspan_t facet_normals(facet_normalsb.data(), n_shape); // Working memory for loop std::vector coordinate_dofsb(num_dofs_g * gdim); - cmdspan2_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, gdim); + mdspan_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, + gdim); std::array Jb; std::array Kb; mdspan2_t J(Jb.data(), gdim, tdim); @@ -532,7 +533,8 @@ dolfinx_contact::Contact::pack_gap(int pair) const const int tdim = topology->dim(); std::vector coordinate_dofsb(num_dofs_g * gdim); - cmdspan2_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, gdim); + mdspan_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, + gdim); std::array coordb; mdspan2_t coord(coordb.data(), 1, gdim); @@ -1021,11 +1023,12 @@ dolfinx_contact::Contact::pack_ny(int pair) const candidate_mesh->topology()->cell_type()); auto [facet_normalsb, n_shape] = basix::cell::facet_outward_normals(cell_type); - cmdspan2_t facet_normals(facet_normalsb.data(), n_shape); + mdspan_t facet_normals(facet_normalsb.data(), n_shape); // Working memory for loop std::vector coordinate_dofsb(num_dofs_g * gdim); - cmdspan2_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, gdim); + mdspan_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, + gdim); std::array Jb; std::array Kb; mdspan2_t J(Jb.data(), gdim, tdim); diff --git a/cpp/KernelData.cpp b/cpp/KernelData.cpp index 85613ece..4fdc0616 100644 --- a/cpp/KernelData.cpp +++ b/cpp/KernelData.cpp @@ -92,7 +92,8 @@ dolfinx_contact::KernelData::KernelData( double dolfinx_contact::KernelData::compute_first_facet_jacobian( const std::size_t facet_index, dolfinx_contact::mdspan2_t J, dolfinx_contact::mdspan2_t K, dolfinx_contact::mdspan2_t J_tot, - std::span detJ_scratch, dolfinx_contact::cmdspan2_t coords) const + std::span detJ_scratch, + dolfinx_contact::mdspan_t coords) const { dolfinx_contact::cmdspan4_t full_basis(_c_basis_values.data(), _c_basis_shape); @@ -112,12 +113,13 @@ double dolfinx_contact::KernelData::compute_first_facet_jacobian( } //----------------------------------------------------------------------------- void dolfinx_contact::KernelData::update_normal( - std::span n, dolfinx_contact::cmdspan2_t K, + std::span n, dolfinx_contact::mdspan_t K, const std::size_t local_index) const { - return _update_normal( - n, K, dolfinx_contact::cmdspan2_t(_facet_normals.data(), _normals_shape), - local_index); + return _update_normal(n, K, + dolfinx_contact::mdspan_t( + _facet_normals.data(), _normals_shape), + local_index); } //----------------------------------------------------------------------------- std::span diff --git a/cpp/KernelData.h b/cpp/KernelData.h index b3bb3331..fdb79c13 100644 --- a/cpp/KernelData.h +++ b/cpp/KernelData.h @@ -20,12 +20,13 @@ namespace dolfinx_contact { -using jac_fn = std::function, cmdspan2_t, s_cmdspan2_t, - cmdspan2_t)>; +using jac_fn = std::function, + mdspan_t, s_cmdspan2_t, mdspan_t)>; -using normal_fn = std::function, cmdspan2_t, cmdspan2_t, - const std::size_t)>; +using normal_fn + = std::function, mdspan_t, + mdspan_t, const std::size_t)>; class KernelData { @@ -100,9 +101,9 @@ class KernelData const std::vector& offsets_array() const { return _offsets; } // Return reference facet normals - cmdspan2_t facet_normals() const + mdspan_t facet_normals() const { - return cmdspan2_t(_facet_normals.data(), _normals_shape); + return mdspan_t(_facet_normals.data(), _normals_shape); } /// Compute the following jacobians on a given facet: @@ -119,7 +120,7 @@ class KernelData double update_jacobian(std::size_t q, std::size_t facet_index, double detJ, mdspan2_t J, mdspan2_t K, mdspan2_t J_tot, std::span detJ_scratch, - cmdspan2_t coords) const + mdspan_t coords) const { cmdspan4_t full_basis(_c_basis_values.data(), _c_basis_shape); const std::size_t q_pos = _qp_offsets[facet_index] + q; @@ -148,13 +149,13 @@ class KernelData double compute_first_facet_jacobian(std::size_t facet_index, mdspan2_t J, mdspan2_t K, mdspan2_t J_tot, std::span detJ_scratch, - cmdspan2_t coords) const; + mdspan_t coords) const; /// update normal /// @param[in, out] n The facet normal /// @param[in] K The inverse Jacobian /// @param[in] local_index The facet index local to the cell - void update_normal(std::span n, cmdspan2_t K, + void update_normal(std::span n, mdspan_t K, std::size_t local_index) const; /// Return quadrature weights for the i-th facet diff --git a/cpp/QuadratureRule.cpp b/cpp/QuadratureRule.cpp index a2ba1299..2af657a4 100644 --- a/cpp/QuadratureRule.cpp +++ b/cpp/QuadratureRule.cpp @@ -34,7 +34,8 @@ dolfinx_contact::QuadratureRule::QuadratureRule(dolfinx::mesh::CellType ct, std::vector& q_weights = quadrature.back(); std::size_t num_points = q_weights.size(); std::size_t pt_shape = quadrature.front().size() / num_points; - cmdspan2_t qp(quadrature.front().data(), num_points, pt_shape); + mdspan_t qp(quadrature.front().data(), num_points, + pt_shape); _points = std::vector(num_points * _num_sub_entities * _tdim); _entity_offset = std::vector(_num_sub_entities + 1, 0); @@ -91,7 +92,7 @@ dolfinx_contact::QuadratureRule::QuadratureRule(dolfinx::mesh::CellType ct, auto [sub_geomb, sub_geom_shape] = basix::cell::sub_entity_geometry(b_ct, dim, i); - cmdspan2_t coords(sub_geomb.data(), sub_geom_shape); + mdspan_t coords(sub_geomb.data(), sub_geom_shape); // Push forward quadrature point from reference entity to reference // entity on cell @@ -132,10 +133,10 @@ std::size_t QuadratureRule::num_points(int i) const return _entity_offset[i + 1] - _entity_offset[i]; } //----------------------------------------------------------------------------------------------- -cmdspan2_t QuadratureRule::points(int i) const +mdspan_t QuadratureRule::points(int i) const { assert(i < _num_sub_entities); - cmdspan2_t all_points(_points.data(), _weights.size(), _tdim); + mdspan_t all_points(_points.data(), _weights.size(), _tdim); return MDSPAN_IMPL_STANDARD_NAMESPACE::submdspan( all_points, std::pair(_entity_offset[i], _entity_offset[i + 1]), MDSPAN_IMPL_STANDARD_NAMESPACE::full_extent); diff --git a/cpp/QuadratureRule.h b/cpp/QuadratureRule.h index 65e8727f..46789972 100644 --- a/cpp/QuadratureRule.h +++ b/cpp/QuadratureRule.h @@ -21,8 +21,6 @@ using mdspan_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< using mdspan2_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< double, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents>; -using cmdspan2_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< - const double, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents>; using mdspan3_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< double, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents>; using mdspan4_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< @@ -35,9 +33,6 @@ using s_cmdspan2_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< using s_cmdspan3_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< const double, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents, stdex::layout_stride>; -using s_cmdspan4_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< - const double, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents, - stdex::layout_stride>; class QuadratureRule { @@ -85,7 +80,7 @@ class QuadratureRule /// Return the quadrature points for the ith entity /// @param[in] i The local entity index - cmdspan2_t points(int i) const; + mdspan_t points(int i) const; /// Return the quadrature weights for the ith entity /// @param[in] i The local entity index diff --git a/cpp/RayTracing.h b/cpp/RayTracing.h index 58380d95..b43b61c6 100644 --- a/cpp/RayTracing.h +++ b/cpp/RayTracing.h @@ -273,7 +273,7 @@ int raytracing_cell( auto dphi = MDSPAN_IMPL_STANDARD_NAMESPACE::submdspan( basis, std::pair{1, tdim + 1}, 0, MDSPAN_IMPL_STANDARD_NAMESPACE::full_extent, 0); - cmdspan2_t coords(coordinate_dofs.data(), cmap.dim(), gdim); + mdspan_t coords(coordinate_dofs.data(), cmap.dim(), gdim); mdspan2_t _xk(x_k.data(), 1, gdim); for (int k = 0; k < max_iter; ++k) { @@ -500,7 +500,7 @@ compute_ray(const dolfinx::mesh::Mesh& mesh, std::span xi, std::span X) { const std::vector& facet = facets[facet_index]; - dolfinx_contact::cmdspan2_t x(xb.data(), x_shape); + dolfinx_contact::mdspan_t x(xb.data(), x_shape); for (std::size_t i = 0; i < tdim; ++i) { X[i] = x(facet.front(), i); diff --git a/cpp/coefficients.cpp b/cpp/coefficients.cpp index 7e97a671..34c96592 100644 --- a/cpp/coefficients.cpp +++ b/cpp/coefficients.cpp @@ -15,9 +15,9 @@ using namespace dolfinx_contact; void dolfinx_contact::transformed_push_forward( const dolfinx::fem::FiniteElement* element, cmdspan4_t reference_basis, std::vector& element_basisb, - mdspan3_t basis_values, cmdspan2_t J, cmdspan2_t K, double detJ, - std::size_t basis_offset, std::size_t q, std::int32_t cell, - std::span cell_info) + mdspan3_t basis_values, mdspan_t J, + mdspan_t K, double detJ, std::size_t basis_offset, + std::size_t q, std::int32_t cell, std::span cell_info) { const std::function&, const std::span&, std::int32_t, @@ -28,8 +28,10 @@ void dolfinx_contact::transformed_push_forward( // Get push forward function auto push_forward_fn = element->basix_element() - .map_fn(); + .map_fn, + dolfinx_contact::mdspan_t, + dolfinx_contact::mdspan_t>(); mdspan2_t element_basis(element_basisb.data(), basis_values.extent(1), basis_values.extent(2)); diff --git a/cpp/coefficients.h b/cpp/coefficients.h index cf0958fb..daa106a1 100644 --- a/cpp/coefficients.h +++ b/cpp/coefficients.h @@ -30,9 +30,9 @@ namespace dolfinx_contact void transformed_push_forward( const dolfinx::fem::FiniteElement* element, cmdspan4_t reference_basis, std::vector& element_basisb, - mdspan3_t basis_values, cmdspan2_t J, cmdspan2_t K, double detJ, - std::size_t basis_offset, std::size_t q, std::int32_t cell, - std::span cell_info); + mdspan3_t basis_values, mdspan_t J, + mdspan_t K, double detJ, std::size_t basis_offset, + std::size_t q, std::int32_t cell, std::span cell_info); /// @brief Pack a coefficient at quadrature points. /// diff --git a/cpp/contact_kernels.cpp b/cpp/contact_kernels.cpp index 26fe9806..c544fcac 100644 --- a/cpp/contact_kernels.cpp +++ b/cpp/contact_kernels.cpp @@ -72,7 +72,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -221,7 +222,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -387,7 +389,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -571,7 +574,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -802,7 +806,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -999,7 +1004,8 @@ dolfinx_contact::generate_contact_kernel( const std::uint32_t tdim = kd.tdim(); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler diff --git a/cpp/elasticity.cpp b/cpp/elasticity.cpp index 3e077b32..5402bab7 100644 --- a/cpp/elasticity.cpp +++ b/cpp/elasticity.cpp @@ -8,7 +8,7 @@ //----------------------------------------------------------------------------- void dolfinx_contact::compute_normal_strain_basis( - mdspan2_t epsn, mdspan2_t tr, dolfinx_contact::cmdspan2_t K, + mdspan2_t epsn, mdspan2_t tr, dolfinx_contact::mdspan_t K, dolfinx_contact::s_cmdspan3_t dphi, std::array n_surf, std::span n_phys, std::size_t q_pos) { @@ -42,11 +42,9 @@ void dolfinx_contact::compute_normal_strain_basis( } } //----------------------------------------------------------------------------- -void dolfinx_contact::compute_sigma_n_basis(mdspan3_t sig_n, cmdspan2_t K, - s_cmdspan3_t dphi, - std::span n, - double mu, double lmbda, - std::size_t q_pos) +void dolfinx_contact::compute_sigma_n_basis( + mdspan3_t sig_n, mdspan_t K, s_cmdspan3_t dphi, + std::span n, double mu, double lmbda, std::size_t q_pos) { const std::size_t ndofs_cell = sig_n.extent(0); const std::size_t bs = K.extent(1); diff --git a/cpp/elasticity.h b/cpp/elasticity.h index 72f213ab..7453e8a9 100644 --- a/cpp/elasticity.h +++ b/cpp/elasticity.h @@ -27,7 +27,8 @@ namespace dolfinx_contact /// @param[in] n_1 1st normal vector, typically n_surf /// @param[in] n_2 2nd normal vector, typically n_phys /// @param[in] q_pos offset of quadrature point for accessing dphi -void compute_normal_strain_basis(mdspan2_t epsn, mdspan2_t tr, cmdspan2_t K, +void compute_normal_strain_basis(mdspan2_t epsn, mdspan2_t tr, + mdspan_t K, s_cmdspan3_t s_dphi, std::array n_1, std::span n_2, std::size_t q_pos); @@ -46,9 +47,9 @@ void compute_normal_strain_basis(mdspan2_t epsn, mdspan2_t tr, cmdspan2_t K, /// @param[in] mu The poisson ratio /// @param[in] lmbda The 1st Lame parameter /// @param[in] q_pos The offset of quadrature point for accessing dphi -void compute_sigma_n_basis(mdspan3_t sig_n, cmdspan2_t K, s_cmdspan3_t dphi, - std::span n, double mu, double lmbda, - std::size_t q_pos); +void compute_sigma_n_basis(mdspan3_t sig_n, mdspan_t K, + s_cmdspan3_t dphi, std::span n, + double mu, double lmbda, std::size_t q_pos); /// @brief Compute sigma(u)*n from grad(u) /// diff --git a/cpp/geometric_quantities.cpp b/cpp/geometric_quantities.cpp index b15477cb..f616cbaf 100644 --- a/cpp/geometric_quantities.cpp +++ b/cpp/geometric_quantities.cpp @@ -21,7 +21,7 @@ void dolfinx_contact::pull_back_nonaffine( std::span X, std::span work_array, std::span x, const dolfinx::fem::CoordinateElement& cmap, - cmdspan2_t cell_geometry, double tol, const int max_it) + mdspan_t cell_geometry, double tol, const int max_it) { assert((std::size_t)cmap.dim() == cell_geometry.extent(0)); // Temporary data structures for Newton iteration @@ -95,9 +95,10 @@ void dolfinx_contact::pull_back_nonaffine( std::array dolfinx_contact::push_forward_facet_normal( std::span work_array, std::span x, std::size_t gdim, - std::size_t tdim, cmdspan2_t coordinate_dofs, const std::size_t facet_index, + std::size_t tdim, mdspan_t coordinate_dofs, + const std::size_t facet_index, const dolfinx::fem::CoordinateElement& cmap, - cmdspan2_t reference_normals) + mdspan_t reference_normals) { const std::array c_shape = cmap.tabulate_shape(1, 1); const std::size_t basis_size @@ -156,7 +157,8 @@ std::array dolfinx_contact::push_forward_facet_normal( //----------------------------------------------------------------------------- double dolfinx_contact::compute_circumradius(const dolfinx::mesh::Mesh& mesh, - double detJ, cmdspan2_t coordinate_dofs) + double detJ, + mdspan_t coordinate_dofs) { const dolfinx::mesh::CellType cell_type = mesh.topology()->cell_type(); const int gdim = mesh.geometry().dim(); diff --git a/cpp/geometric_quantities.h b/cpp/geometric_quantities.h index 0e33f2a9..4545b023 100644 --- a/cpp/geometric_quantities.h +++ b/cpp/geometric_quantities.h @@ -31,8 +31,10 @@ namespace dolfinx_contact /// @param[in] reference_normals: The facet normals on the reference cell std::array push_forward_facet_normal( std::span work_array, std::span x, std::size_t gdim, - std::size_t tdim, cmdspan2_t coordinate_dofs, const std::size_t facet_index, - const dolfinx::fem::CoordinateElement& cmap, cmdspan2_t reference_normals); + std::size_t tdim, mdspan_t coordinate_dofs, + const std::size_t facet_index, + const dolfinx::fem::CoordinateElement& cmap, + mdspan_t reference_normals); /// @brief Allocate memory for pull-back on a non affine cell for a /// single point. @@ -40,9 +42,8 @@ std::array push_forward_facet_normal( /// @param[in] gdim The geometrical dimension /// @param[in] tdim The topological dimension /// @returns Vector of sufficient size -std::vector -allocate_pull_back_nonaffine(const dolfinx::fem::CoordinateElement& cmap, - int gdim, int tdim); +std::vector allocate_pull_back_nonaffine( + const dolfinx::fem::CoordinateElement& cmap, int gdim, int tdim); /// @brief Pull back a single point of a non-affine cell to the /// reference cell. @@ -61,8 +62,8 @@ allocate_pull_back_nonaffine(const dolfinx::fem::CoordinateElement& cmap void pull_back_nonaffine(std::span X, std::span work_array, std::span x, const dolfinx::fem::CoordinateElement& cmap, - cmdspan2_t cell_geometry, double tol = 1e-8, - const int max_it = 10); + mdspan_t cell_geometry, + double tol = 1e-8, const int max_it = 10); /// Compute circumradius for a cell with given coordinates and /// determinant of Jacobian @@ -72,7 +73,8 @@ void pull_back_nonaffine(std::span X, std::span work_array, /// @param[in] coordinate_dofs The cell geometry /// @returns The circumradius of the cell double compute_circumradius(const dolfinx::mesh::Mesh& mesh, - double detJ, cmdspan2_t coordinate_dofs); + double detJ, + mdspan_t coordinate_dofs); /// @brief Push forward facet normal /// diff --git a/cpp/meshtie_kernels.cpp b/cpp/meshtie_kernels.cpp index 9e1cc650..29be0b3b 100644 --- a/cpp/meshtie_kernels.cpp +++ b/cpp/meshtie_kernels.cpp @@ -44,7 +44,8 @@ dolfinx_contact::generate_meshtie_kernel( auto tdim = std::size_t(kd.tdim()); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -186,7 +187,8 @@ dolfinx_contact::generate_meshtie_kernel( auto tdim = std::size_t(kd.tdim()); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -286,7 +288,8 @@ dolfinx_contact::generate_meshtie_kernel( // Retrieve some data from kd auto tdim = std::size_t(kd.tdim()); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -483,7 +486,8 @@ dolfinx_contact::generate_poisson_kernel( auto tdim = std::size_t(kd.tdim()); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler @@ -609,7 +613,8 @@ dolfinx_contact::generate_poisson_kernel( // Retrieve some data from kd auto tdim = std::size_t(kd.tdim()); // NOTE: DOLFINx has 3D input coordinate dofs - cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), 3); + mdspan_t coord(coordinate_dofs, kd.num_coordinate_dofs(), + 3); // Create data structures for jacobians // We allocate more memory than required, but its better for the compiler diff --git a/cpp/rigid_surface_kernels.cpp b/cpp/rigid_surface_kernels.cpp index aa9f89fb..c6f8b22d 100644 --- a/cpp/rigid_surface_kernels.cpp +++ b/cpp/rigid_surface_kernels.cpp @@ -75,8 +75,8 @@ dolfinx_contact::generate_rigid_surface_kernel( const std::size_t ndofs_cell = kd.ndofs_cell(); // Reshape coordinate dofs to two-dimensional array - dolfinx_contact::cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), - 3); + dolfinx_contact::mdspan_t coord( + coordinate_dofs, kd.num_coordinate_dofs(), 3); // Compute Jacobian and determinant at first quadrature point std::array Jb; @@ -234,8 +234,8 @@ dolfinx_contact::generate_rigid_surface_kernel( const std::uint32_t ndofs_cell = kd.ndofs_cell(); // Reshape coordinate dofs to two-dimensional array - dolfinx_contact::cmdspan2_t coord(coordinate_dofs, kd.num_coordinate_dofs(), - 3); + dolfinx_contact::mdspan_t coord( + coordinate_dofs, kd.num_coordinate_dofs(), 3); // Compute Jacobian and determinant at first quadrature point std::array Jb; diff --git a/cpp/utils.cpp b/cpp/utils.cpp index a289046f..cc215c70 100644 --- a/cpp/utils.cpp +++ b/cpp/utils.cpp @@ -65,8 +65,9 @@ dolfinx_contact::read_mesh(const std::string& filename, //----------------------------------------------------------------------------- void dolfinx_contact::pull_back( dolfinx_contact::mdspan3_t J, dolfinx_contact::mdspan3_t K, - std::span detJ, std::span X, dolfinx_contact::cmdspan2_t x, - dolfinx_contact::cmdspan2_t coordinate_dofs, + std::span detJ, std::span X, + dolfinx_contact::mdspan_t x, + dolfinx_contact::mdspan_t coordinate_dofs, const dolfinx::fem::CoordinateElement& cmap) { const std::size_t num_points = x.extent(0); @@ -566,8 +567,9 @@ void dolfinx_contact::evaluate_basis_functions( double dolfinx_contact::compute_facet_jacobian( dolfinx_contact::mdspan2_t J, dolfinx_contact::mdspan2_t K, dolfinx_contact::mdspan2_t J_tot, std::span detJ_scratch, - dolfinx_contact::cmdspan2_t J_f, dolfinx_contact::s_cmdspan2_t dphi, - dolfinx_contact::cmdspan2_t coords) + dolfinx_contact::mdspan_t J_f, + dolfinx_contact::s_cmdspan2_t dphi, + dolfinx_contact::mdspan_t coords) { std::size_t gdim = J.extent(0); auto coordinate_dofs = MDSPAN_IMPL_STANDARD_NAMESPACE::submdspan( @@ -589,33 +591,36 @@ double dolfinx_contact::compute_facet_jacobian( //------------------------------------------------------------------------------------- std::function, dolfinx_contact::cmdspan2_t, - dolfinx_contact::s_cmdspan2_t, dolfinx_contact::cmdspan2_t)> + dolfinx_contact::mdspan2_t, std::span, + dolfinx_contact::mdspan_t, dolfinx_contact::s_cmdspan2_t, + dolfinx_contact::mdspan_t)> dolfinx_contact::get_update_jacobian_dependencies( const dolfinx::fem::CoordinateElement& cmap) { if (cmap.is_affine()) { // Return function that returns the input determinant - return [](double detJ, [[maybe_unused]] dolfinx_contact::mdspan2_t J, - [[maybe_unused]] dolfinx_contact::mdspan2_t K, - [[maybe_unused]] dolfinx_contact::mdspan2_t J_tot, - [[maybe_unused]] std::span detJ_scratch, - [[maybe_unused]] dolfinx_contact::cmdspan2_t J_f, - [[maybe_unused]] dolfinx_contact::s_cmdspan2_t dphi, - [[maybe_unused]] dolfinx_contact::cmdspan2_t coords) + return + [](double detJ, [[maybe_unused]] dolfinx_contact::mdspan2_t J, + [[maybe_unused]] dolfinx_contact::mdspan2_t K, + [[maybe_unused]] dolfinx_contact::mdspan2_t J_tot, + [[maybe_unused]] std::span detJ_scratch, + [[maybe_unused]] dolfinx_contact::mdspan_t J_f, + [[maybe_unused]] dolfinx_contact::s_cmdspan2_t dphi, + [[maybe_unused]] dolfinx_contact::mdspan_t coords) { return detJ; }; } else { // Return function that returns the input determinant - return [](double detJ, [[maybe_unused]] dolfinx_contact::mdspan2_t J, - [[maybe_unused]] dolfinx_contact::mdspan2_t K, - [[maybe_unused]] dolfinx_contact::mdspan2_t J_tot, - [[maybe_unused]] std::span detJ_scratch, - [[maybe_unused]] dolfinx_contact::cmdspan2_t J_f, - [[maybe_unused]] dolfinx_contact::s_cmdspan2_t dphi, - [[maybe_unused]] dolfinx_contact::cmdspan2_t coords) + return + [](double detJ, [[maybe_unused]] dolfinx_contact::mdspan2_t J, + [[maybe_unused]] dolfinx_contact::mdspan2_t K, + [[maybe_unused]] dolfinx_contact::mdspan2_t J_tot, + [[maybe_unused]] std::span detJ_scratch, + [[maybe_unused]] dolfinx_contact::mdspan_t J_f, + [[maybe_unused]] dolfinx_contact::s_cmdspan2_t dphi, + [[maybe_unused]] dolfinx_contact::mdspan_t coords) { double new_detJ = dolfinx_contact::compute_facet_jacobian( J, K, J_tot, detJ_scratch, J_f, dphi, coords); @@ -624,8 +629,9 @@ dolfinx_contact::get_update_jacobian_dependencies( } } //------------------------------------------------------------------------------------- -std::function, dolfinx_contact::cmdspan2_t, - dolfinx_contact::cmdspan2_t, const std::size_t)> +std::function< + void(std::span, dolfinx_contact::mdspan_t, + dolfinx_contact::mdspan_t, const std::size_t)> dolfinx_contact::get_update_normal( const dolfinx::fem::CoordinateElement& cmap) { @@ -633,8 +639,8 @@ dolfinx_contact::get_update_normal( { // Return function that returns the input determinant return []([[maybe_unused]] std::span n, - [[maybe_unused]] dolfinx_contact::cmdspan2_t K, - [[maybe_unused]] dolfinx_contact::cmdspan2_t n_ref, + [[maybe_unused]] dolfinx_contact::mdspan_t K, + [[maybe_unused]] dolfinx_contact::mdspan_t n_ref, [[maybe_unused]] const std::size_t local_index) { // Do nothing @@ -643,8 +649,9 @@ dolfinx_contact::get_update_normal( else { // Return function that updates the physical normal based on K - return [](std::span n, dolfinx_contact::cmdspan2_t K, - dolfinx_contact::cmdspan2_t n_ref, const std::size_t local_index) + return [](std::span n, dolfinx_contact::mdspan_t K, + dolfinx_contact::mdspan_t n_ref, + const std::size_t local_index) { std::fill(n.begin(), n.end(), 0); auto n_f = MDSPAN_IMPL_STANDARD_NAMESPACE::submdspan( @@ -948,8 +955,8 @@ void dolfinx_contact::compute_physical_points( // Temporary data array std::vector coordinate_dofsb(num_dofs_g * gdim); - dolfinx_contact::cmdspan2_t coordinate_dofs(coordinate_dofsb.data(), - num_dofs_g, gdim); + dolfinx_contact::mdspan_t coordinate_dofs( + coordinate_dofsb.data(), num_dofs_g, gdim); for (std::size_t i = 0; i < facets.size(); i += 2) { auto x_dofs = MDSPAN_IMPL_STANDARD_NAMESPACE::submdspan( diff --git a/cpp/utils.h b/cpp/utils.h index e499cb58..4e7658b1 100644 --- a/cpp/utils.h +++ b/cpp/utils.h @@ -98,7 +98,8 @@ using kernel_fn /// @param[in] coordinate_dofs: geometry coordinates of cell /// @param[in] cmap: the coordinate element void pull_back(mdspan3_t J, mdspan3_t K, std::span detJ, - std::span X, cmdspan2_t x, cmdspan2_t coordinate_dofs, + std::span X, mdspan_t x, + mdspan_t coordinate_dofs, const dolfinx::fem::CoordinateElement& cmap); /// @param[in] cells: the cells to be sorted @@ -211,8 +212,9 @@ void evaluate_basis_functions(const dolfinx::fem::FunctionSpace& V, /// @param[in] coords - the coordinates of the facet /// @return absolute value of determinant of J_tot double compute_facet_jacobian(mdspan2_t J, mdspan2_t K, mdspan2_t J_tot, - std::span detJ_scratch, cmdspan2_t J_f, - s_cmdspan2_t dphi, cmdspan2_t coords); + std::span detJ_scratch, + mdspan_t J_f, s_cmdspan2_t dphi, + mdspan_t coords); /// @brief Convenience function to update Jacobians /// @@ -221,7 +223,8 @@ double compute_facet_jacobian(mdspan2_t J, mdspan2_t K, mdspan2_t J_tot, /// (J*J_f) is computed. /// @param[in] cmap The coordinate element std::function, - cmdspan2_t, s_cmdspan2_t, cmdspan2_t)> + mdspan_t, s_cmdspan2_t, + mdspan_t)> get_update_jacobian_dependencies( const dolfinx::fem::CoordinateElement& cmap); @@ -231,8 +234,8 @@ get_update_jacobian_dependencies( /// For non-affine geometries, a function updating the physical facet normal is /// returned. /// @param[in] cmap The coordinate element -std::function, cmdspan2_t, cmdspan2_t, - const std::size_t)> +std::function, mdspan_t, + mdspan_t, const std::size_t)> get_update_normal(const dolfinx::fem::CoordinateElement& cmap); /// @brief Convert local entity indices to integration entities @@ -478,7 +481,8 @@ compute_projection_map(const dolfinx::mesh::Mesh& mesh, MDSPAN_IMPL_STANDARD_NAMESPACE::dextents> x_dofmap = mesh.geometry().dofmap(); std::vector coordinate_dofsb(num_dofs_g * gdim); - cmdspan2_t coordinate_dofs(coordinate_dofsb.data(), num_dofs_g, gdim); + mdspan_t coordinate_dofs(coordinate_dofsb.data(), + num_dofs_g, gdim); auto f_to_c = mesh.topology()->connectivity(tdim - 1, tdim); if (!f_to_c) throw std::runtime_error("Missing facet to cell connectivity"); @@ -506,8 +510,8 @@ compute_projection_map(const dolfinx::mesh::Mesh& mesh, // in a single cell at the same time // Pull back coordinates std::fill(Jb.begin(), Jb.end(), 0); - pull_back(J, K, detJ, X, cmdspan2_t(x.data(), 1, gdim), coordinate_dofs, - cmap); + pull_back(J, K, detJ, X, mdspan_t(x.data(), 1, gdim), + coordinate_dofs, cmap); // Copy into output std::copy_n(X.begin(), tdim, std::next(candidate_X.begin(), i * tdim)); } @@ -593,7 +597,8 @@ compute_raytracing_map(const dolfinx::mesh::Mesh& quadrature_mesh, q_dofmap = geom_q.dofmap(); const std::size_t num_nodes_q = cmap_q.dim(); std::vector coordinate_dofs_qb(num_nodes_q * gdim); - cmdspan2_t coordinate_dofs_q(coordinate_dofs_qb.data(), num_nodes_q, gdim); + mdspan_t coordinate_dofs_q(coordinate_dofs_qb.data(), + num_nodes_q, gdim); auto [reference_normals, rn_shape] = basix::cell::facet_outward_normals( dolfinx::mesh::cell_type_to_basix_type(top_q->cell_type())); @@ -743,7 +748,7 @@ compute_raytracing_map(const dolfinx::mesh::Mesh& quadrature_mesh, std::span X) { const std::vector& facet = bfacets[facet_index_c]; - dolfinx_contact::cmdspan2_t x(xb.data(), x_shape); + dolfinx_contact::mdspan_t x(xb.data(), x_shape); const int f0 = facet.front(); for (std::size_t i = 0; i < tdim; ++i) { diff --git a/python/dolfinx_contact/wrappers.cpp b/python/dolfinx_contact/wrappers.cpp index 1997efa4..6032af4b 100644 --- a/python/dolfinx_contact/wrappers.cpp +++ b/python/dolfinx_contact/wrappers.cpp @@ -87,7 +87,8 @@ NB_MODULE(cpp, m) .def("points", [](dolfinx_contact::QuadratureRule& self, int i) { - dolfinx_contact::cmdspan2_t points_i = self.points(i); + dolfinx_contact::mdspan_t points_i + = self.points(i); std::array shape = {points_i.extent(0), points_i.extent(1)}; std::vector _points(shape[0] * shape[1]);