| 228 | /// @return A mesh geometry. |
| 229 | template <typename U> |
| 230 | Geometry<typename std::remove_reference_t<typename U::value_type>> |
| 231 | create_geometry(const Topology& topology, |
| 232 | const std::vector<fem::CoordinateElement< |
| 233 | std::remove_reference_t<typename U::value_type>>>& elements, |
| 234 | std::span<const std::int64_t> nodes, |
| 235 | std::span<const std::int64_t> xdofs, const U& x, int dim, |
| 236 | const std::function<std::vector<int>( |
| 237 | const graph::AdjacencyList<std::int32_t>&)>& reorder_fn |
| 238 | = nullptr) |
| 239 | { |
| 240 | spdlog::info("Create Geometry (multiple)"); |
| 241 | |
| 242 | assert(std::ranges::is_sorted(nodes)); |
| 243 | using T = typename std::remove_reference_t<typename U::value_type>; |
| 244 | |
| 245 | // Check elements match cell types in topology |
| 246 | const int tdim = topology.dim(); |
| 247 | const std::size_t num_cell_types = topology.entity_types(tdim).size(); |
| 248 | if (elements.size() != num_cell_types) |
| 249 | throw std::runtime_error("Mismatch between topology and geometry."); |
| 250 | |
| 251 | std::vector<fem::ElementDofLayout> dof_layouts; |
| 252 | dof_layouts.reserve(elements.size()); |
| 253 | for (auto& el : elements) |
| 254 | dof_layouts.push_back(el.create_dof_layout()); |
| 255 | |
| 256 | spdlog::info("Got {} dof layouts", dof_layouts.size()); |
| 257 | |
| 258 | // Build 'geometry' dofmap on the topology |
| 259 | auto [_dof_index_map, bs, dofmaps] |
| 260 | = fem::build_dofmap_data(topology.index_maps(topology.dim())[0]->comm(), |
| 261 | topology, dof_layouts, reorder_fn); |
| 262 | auto dof_index_map |
| 263 | = std::make_shared<common::IndexMap>(std::move(_dof_index_map)); |
| 264 | |
| 265 | // If the mesh has higher order geometry, permute the dofmap |
| 266 | if (elements.front().needs_dof_permutations()) |
| 267 | { |
| 268 | const std::int32_t num_cells |
| 269 | = topology.connectivity(topology.dim(), 0)->num_nodes(); |
| 270 | const std::vector<std::uint32_t>& cell_info |
| 271 | = topology.get_cell_permutation_info(); |
| 272 | int d = elements.front().dim(); |
| 273 | for (std::int32_t cell = 0; cell < num_cells; ++cell) |
| 274 | { |
| 275 | std::span dofs(dofmaps.front().data() + cell * d, d); |
| 276 | elements.front().permute_inv(dofs, cell_info[cell]); |
| 277 | } |
| 278 | } |
| 279 | |
| 280 | spdlog::info("Calling compute_local_to_global"); |
| 281 | // Compute local-to-global map from local indices in dofmap to the |
| 282 | // corresponding global indices in cells, and pass to function to |
| 283 | // compute local (dof) to local (position in coords) map from (i) |
| 284 | // local-to-global for dofs and (ii) local-to-global for entries in |
| 285 | // coords |
| 286 | |
| 287 | spdlog::info("xdofs.size = {}", xdofs.size()); |
no test coverage detected